New Spectrophotometric Method for Quantitative Characterization of Density-Driven Convective Instability

Abstract

CO₂ convective dissolution has been regarded as one of the fundamental mechanisms to accelerate the mass transfer of CO₂ into brine. We present a new spectrophotometric method to characterize the convective instability and measure the dissolved CO₂ mass, which enables the real-time quantitative visualization of CO₂/brine transport mechanisms. Successive images were captured to identify the finger development regimes, and the convection morphologies were analyzed by the fingers length and affected area. CO₂ solubility was experimentally studied, and the results are in agreement with the theoretical calculations. CO₂ mass transfer flux was investigated as the Sherwood number changed. The increase in salinity and temperature has a negative effect on CO₂ dissolution; here, numerical simulation and experimental phenomena are qualitatively consistent. In general, these findings confirm the feasibility of the method and improve the understanding of the physical process of CO₂ convective dissolution, which can help assess the CO₂ solubility trapping mass.

1. Introduction

The Earth’s climate changes as a result of the increase of greenhouse gas emissions, particularly emissions of carbon dioxide (CO₂) from burning fossil fuels [1]. Carbon capture and storage (CCS) technology is often considered as a cost-effective solution to mitigate climate change [2]. Saline aquifers have been recognized as promising storage sites due to their large capacity and widespread distribution across the world [3]. One fundamental physical phenomenon related to CO₂ geological storage in saline aquifers is CO₂ dissolution into brine [4]. Once CO₂ is injected in brine-bearing formations, the less dense CO₂ rises upwards until it is confined by an impermeable caprock where it spreads laterally beneath; then, it starts to dissolve into the underlying brine [5]. CO₂-dissolved brine is 0.1–1% denser than the original brine [6]. Thus, a gravitationally unstable system develops with CO₂-dissolved brine overlying the original brine, and this leads to density-driven convection and triggers instabilities.

CO₂ dissolution in brine causes density-driven convection, which is a key mechanism for the efficient trapping of CO₂. CO₂ convective dissolution accelerates the CO₂ mass transfer rate and reduces the risk of leakage, which is favorable for long-term sequestration [7]. Although there are some CCS pilot projects distributed in different parts of the world, there is still a lack of relevant knowledge concerning convective dissolution in CCS projects. Furthermore, understanding the convective dissolution process to predict the long-term field evolution of CO₂ is important to sequestration technology implementation [8].

During the CO₂ storage in saline aquifers, the buoyancy of the original formation water is controlled by both thermal (temperature) and compositional (salinity) effects [9]. It is because the solubility of CO₂ is sensitive to the temperature and salinity conditions of saline aquifers. In those actual ongoing CCS projects, the temperature and salinity conditions vary over a wide range [10]. “Cold” sedimentary basins of cold regions and offshore, with the temperature during brine—CO₂ interaction below 50 °C are more favorable for CO₂ storage [11], because CO₂ attains higher density at shallower depths than in “warm” sedimentary basins, which are characterized by high-temperature gradients where dense-fluid conditions are reached at greater depths. In shallow groundwater systems, the salinity of formation water increases with the formation depth, while the interfacial tension (IFT) increases with the increase in salinity, which is attributed to the hydration of ions [12]. Numerical studies have shown that for the constant pressure and temperature, high-salinity brines decrease the density difference between CO₂-dissolved brine and original brine, delaying the onset of instability and slowing down the long-term convective mixing [13].

The instability behavior of convective dissolution has been investigated by many researchers. Lindeberg and Wessel-Berg [14] conducted numerical simulation of the stability criteria for the density-driven convection process. It was followed by several other studies that applied linear stability analysis to predict the convection starting time, the dominant wavelength, and wavenumber for the convective fingers [15,16,17].

As the spatial and temporal scales involved make it difficult to analyze these dynamics in situ, the experimental research on convective dissolution at the laboratory scale is receiving more and more attention [18,19,20]. Three main types of experimental methods have been used to analyze CO₂ convective dissolution in aqueous solutions: PVT (pressure-volume-temperature) reactor experiments [21,22,23], Hele–Shaw cell experiments [7,16,24,25], and bead pack experiments [18,19,26,27]. PVT reactor experimental results can be used to easily quantify the total dissolved CO₂ mass and mass transfer rate at high pressures and elevated temperatures; however, they have limited ability to visualize the convective dissolution process. In contrast, the Hele–Shaw cell experiment is optimal for visual studies, which are commonly used at ambient conditions to visualize the dissolution process; however, the quantification of the mass transfer rate is not as straightforward as that in PVT reactor experiments. In addition, 3D quantitative measurements can be realized by bead pack experiments, but these experiments must be combined with three-dimensional imaging techniques, such as X-ray computed tomography or a magnetic resonance imaging system. In other words, this method needs equipment with complex instruments, and it is not easy to operate.

Previous experiments in Hele–Shaw cells have stimulated considerable theoretical investigations of convective flows. A Hele–Shaw cell is made of two transparent glass or Plexiglas^® plates that are parallel to each other, sealed on the edges, and separated by a narrow gap, which is filled with fluids. Many studies have been conducted with Hele–Shaw cells to reproduce the convective behavior of CO₂-dissolved in brine with analogue fluids [28,29,30,31]. The most common fluids are propylene glycol (PPG)/water, mixtures of methanol and ethylene glycol (EG-MeOH)/water, and potassium permanganate/water. The main purpose of using analogue fluids is to mimic the density differences of the actual fluids and allow experiments to proceed at normal laboratory conditions, length, and time scales [32]. However, there are some fundamental differences between analogue fluids and CO₂/brine systems, such as the density-concentration behavior, viscosity differences, and miscible properties [15,33]. These experimental observations of analogue fluids can aid in understanding the dissolution process and provide physical evidence of mixing induced by density-driven convection in CO₂ geological storage.

Due to the necessity of quantitatively measuring the CO₂ mass transfer to predict CO₂ sequestration capabilities, recent visualization experiments have depended on the optical technique applied using different approaches to quantify the mass transfer. Backhaus et al. measured the temporal evolution using an optical shadowgraph and determined the mass transfer rate [31]. Tsai et al. characterized the maximum convective mass flux by measuring the constant speed of the moving interface [34]. Faisal et al. determined the total mass of dissolved CO₂ by a total carbon analyzer (TC analyzer), which utilizes a catalytic oxidation combustion technique [25]. Rasmusson et al. proposed an experimental method that uses the refractive index of the fluid to quantify the mass flux [35]. However, the previously mentioned methods are either conducted with analogue fluids in CO₂ convective dissolution experiments or cannot measure the CO₂ dissolution mass in a timely and successive manner. There is a limited number of quantitative visual experiments in the area of CO₂ convective dissolution due to the difficulty of measuring the mass transfer quantitatively within the CO₂/brine system. Furthermore, none of the existing studies have compared the effects of temperature and salinity on convective instability through visualization experiments and simulation.

In this study, we carried out in situ dissolved CO₂ concentration measurement experiments using a new spectrophotometric method. We also present a numerical analysis of the density-driven convective instability and the dissolved CO₂ mass in Hele–Shaw cells with respect to different salinity and temperature conditions. The concentration of dissolved CO₂ was reported, and the convective fingers images during the CO₂ dissolution process were quantitatively analyzed. The purpose of this study is to provide quantitative visual evidence of saline water convection with dissolved CO₂ on time and space scales that are easily accessible in the laboratory.

2. Experimental

Experiments were carried out under three different temperature and four salinity conditions. For the cases selected here, temperature conditions can be compared to a storage scenario in a “cold” basin of saline aquifers (geothermal gradient equal to 25 °C/km and average surface temperature around 10 °C) and injection into a rather “shallow” reservoir [36]. Furthermore, sodium chloride (NaCl) was selected to represent the salinity in this study, because Na⁺ is one of the most commonly cations that are found in groundwater [37]. The highest salinity condition, 15,000 mg/L, was inspired by data referring to the CO₂ storage reservoir at Otway Basin [38].

2.1. Experimental Setup

Figure 1 shows a schematic of the experimental setup in this study. The experimental setup consists of a Hele–Shaw cell, syringe pump, CO₂ tank, temperature controller, imaging system, and a PC with a data acquisition system. The two main parts of this setup are the temperature control system and imaging system. It should be pointed out that all experiments were performed at atmospheric pressure.

The glass sheets of the Hele–Shaw cell were separated by spacers with a width b of 1 mm on all sides except at the top. The vertical and bottom sides of the cell were sealed, and the top of the cell was covered but not completely sealed. The length and height of the internal cell dimensions are 100 and 200 mm, respectively. To conduct this experiment, the bottom space of the cell was filled with brine, and the interface height was 100 mm. CO₂ gas was introduced from the upper end, forming the stratified Hele–Shaw flows. Under atmospheric conditions, CO₂ is approximately 1.5 times heavier than air (the molecular weights of CO₂ and air are 44 g/mol and 29 g/mol, respectively); therefore, the air and the excess CO₂ can flow outside of the unsealed top.

The main purpose of Hele–Shaw cells is to simulate an environment that is mathematically analogous to 2D flow in porous media. The flow behavior in the Hele–Shaw cell is governed by the same Darcy’s law as flow in porous media. The permeability of the Hele–Shaw cell was calculated using the equation k = d²/12, and the porosity of the cell was considered as 1.

A custom electro-thermal incubator was used to adjust the temperature of the experiments. It was made with thermal insulation material and an aluminum alloy supported frame; a silicone heating plate was attached to the inner surface. The heating plate related to a power source through the conductive electrode, which was connected with a temperature-control adjuster (AT72AAS4R, Autolise Co., Ltd., Dalian, China, 20–80 °C, resolution 0.1 °C). An ungrounded K-type thermo couple was used to monitor the experimental temperature. The temperature control circuit can provide rapid temperature increases and can maintain a constant predetermined temperature. The front side of the incubator was removable, which was convenient for experimental preparation.

The convective fingers images were captured at regular intervals using an MV-EM510 image processing charge-coupled device (CCD) camera (Microvision Digital Image Technology Co., Ltd., Beijing, China.) with a maximum resolution of 2456 × 2058 pixels and an 8-bit output. This high-resolution camera has a 2/3 progressive scan CCD sensor and can capture 15 frames per second. It was directly connected to a computer to facilitate advanced image processing. A CCD camera was positioned before the Hele–Shaw cell, while a monochromatic light filter and a uniform table light were positioned behind the cell. The central wavelength of the filter is 615 nm, which corresponds to the maximum wavelength of the pH indicator. The CCD relative spectral response near this wavelength is up to 80%, which meets the experimental requirements. Due to the excellent modulation characteristics of the table light, the infused current of the light source was modulated with a light controller to make the incident light brightness change with a signal from the modulator. The light controller is equipped with an RS-232 interface, which can be used for reading and setting light source parameters on the controller panel. The CCD camera, Hele–Shaw cell, light filter, and table light were placed in a row at a proper distance and height so that the camera could focus on the target area. These positions remained the same for subsequent experiments. To avoid external reflections, the above four devices were placed in the custom incubator as a darkroom where the uniform table light was the only source of light.

The pump used for CO₂ injection is an ISCO pump (260D, Teledyne ISCO, Lincoln, NE, USA), which is a precise and accurate laboratory pump with a minimum injection rate of 0.001 mL/min, maximum injection rate of 107 mL/min, and accuracy of 0.5%. This ISCO pump was connected to the CO₂ tank via a tube and a pressure regulator. Before the experiment, it was filled with CO₂ and maintained at atmospheric pressure. The temperature control jacket encircles the cylinder of the pump, allowing water to circulate through the jacket to maintain the temperature of CO₂ within the cylinder.

2.2. Material and Properties

A pH indicator, bromocresol green, was added to the brine for visualizing the dissolved CO₂ convective fingers. Under our experimental conditions, dissolved CO₂ into brine at atmospheric pressure and the solution pH value decreased to approximately 4. The transition range of bromocresol green is 3.8–5.4; therefore, the pH range of the CO₂-dissolved brine was within the pH range of the indicator. To match the intensity of the incident light and obtain the highest quality images, several trial and error experiments were performed to determine the bromocresol green concentration. For our experimental equipment, the optimum concentration of bromocresol green is 2.5 × 10⁻⁴ mol/L. The small amount of bromocresol green has no effect on the properties of the water.

The working solutions were prepared by successive dilution. The pH indicator solution was prepared by dissolving 0.014 g of bromocresol green in 2 mL of 0.025 mol/L sodium hydroxide and adding de-ionized water to obtain a total volume of 200 mL. Then, sodium chloride (NaCl) was dissolved in the pH indicator solution to obtain different concentrations of brine. CO₂ with 99.99% purity was used for experiments.

The density of aqueous CO₂ solutions was calculated according to Duan [39]. The viscosities of the brine solutions were determined by using a rotational viscometer (NDJ-5S Jingtian Co., Ltd., Shanghai, China). The diffusion coefficient was calculated by the Unver and Himmelblau [40] correlation, which related the diffusion coefficients to temperature only.

Table 1. Properties of experimental fluids.

No.	Temperature (°C)	Salinity (mg/L)	∆ρ (kg/m3)	Viscosity (kg/m s)	Diffusion Coefficient (m2/s)	Rayleigh Number
1	25	0	0.37	9.13 × 10−4	1.85 × 10−9	17,406.23
2	25	5000	0.36	9.27 × 10−4	1.85 × 10−9	17,143.36
3	25	10,000	0.35	9.42 × 10−4	1.85 × 10−9	16,401.75
4	25	15,000	0.34	1.03 × 10−3	1.85 × 10−9	15,362.34
5	35	0	0.25	7.08 × 10−4	2.18 × 10−9	13,228.03
6	45	0	0.19	5.94 × 10−4	3.03 × 10−9	8621.23

presents the thermodynamic properties of the fluids used in the experiments.

2.3. Experimental Procedure

The following procedure was used to prepare and observe the density-driven convection of dissolved CO₂ into brine. (1) Prior to each experiment, the Hele–Shaw cell was thoroughly cleaned, and the pH indicator solution and brine were prepared. (2) The cell was partially filled to a specified height (10 cm) with brine containing the pH indicator. Next, the CCD camera, Hele–Shaw cell, light filter, and table light were put into a specific position within the electro-thermal incubator. (3) The Hele–Shaw cells were connected to the ISCO pump via a tube, and a specific incident luminous intensity was obtained by opening and adjusting the infused current of the light. (4) The electro-thermal incubator was closed, and the electro-thermal incubator and water bath were set to the desired temperature and left for approximately half an hour to ensure that the temperature had stabilized. (5) The CCD camera was turned on to capture the dynamics at successive times (0.1 fps), and the valve of the pump was opened to introduce CO₂ into the cell at a constant flow rate (0.5 mL/min). (6) The experiments were repeated three times and continued until the gray value of the images remained constant, which indicated that mixing was complete and the system had equilibrated.

3. Theory and Methods

3.1. Dimensionless Parameters

To better understand the density-driven convection due to CO₂ dissolution, some dimensionless parameters were used to describe and compare the experimental results for different scales and conditions.

The Rayleigh number is the main dimensionless number to describe the gravitational instabilities; it is the ratio between the buoyancy and diffusion forces and was defined as follows:

$$\mathit{Ra} = \frac{\mathit{kgh}\Delta \rho }{\varphi \mu D}$$

(1)

where k is the permeability of the porous medium, $\Delta \rho$ refers to the maximum density difference between the two miscible fluids, h is the characteristic length and equals the height of the interface, μ is the viscosity of the solution, D is the diffusion coefficient, and 𝜙 is the porosity.

The dimensionless flux, characterized by the Sherwood number, Sh, is the ratio of the convection mass flux to the diffusive flux

$$\mathit{Sh} = \frac{F_c}{\varphi \Delta cD/h}$$

(2)

where the convection flux $F_c=\varphi u\Delta c$ is based on the concentration difference between the liquids, $\Delta c$ is the concentration difference, and u is the characteristic velocity.

3.2. Spectrophotometric Method

Interventionary studies involving animals or humans, and other studies that require ethical approval, must list the authority that provided approval and the corresponding ethical approval code.

In this study, we carried out in situ dissolved CO₂ concentration measurement experiments using a spectrophotometric method. According to the Beer–Lambert law [41], light attenuation is related to the properties of the material through which it travels. CO₂ dissolution into saline water was observed under transmitted monochromatic light, and a CCD camera was used to capture images. The gray value of these images are processed using the open-source software Image J [42]; each image contains gray values between 0 (dark) and 255 (bright) and was analyzed to quantify the CO₂ convective dissolution over time.

When CO₂ dissolves in saline water, it reacts with H₂O to form carbonic acid (H₂CO₃), which dissociates into bicarbonate ${\mathrm{HCO}}_3^{-}$ and carbonate ${\mathrm{CO}}_3^{2-}$ ions instantaneously, according to the following equilibrium equations:

$${\mathrm{CO}}_2\left(\mathrm{g}\right) \rightleftharpoons { \mathrm{CO}}_2\left(\mathrm{aq}\right)$$

(3)

$${\mathrm{CO}}_2\left(\mathrm{aq}\right){ + \mathrm{H}}_2\mathrm{O} \rightleftharpoons { \mathrm{H}}^{+}{+ \mathrm{HCO}}_3^{-}$$

(4)

$${\mathrm{HCO}}_3^{-} \rightleftharpoons { \mathrm{H}}^{+}{+ \mathrm{CO}}_3^{2-}.$$

(5)

The equilibrium constant for Equation (4) is several orders of magnitude greater than that for Equation (5). For the convenience of calculation, Equation (4) is neglected. To determine the dissolved CO₂ concentration, a suitable mathematical analysis is used.

Absorbance (A) is the natural logarithm of the ratio of incident to transmitted radiant power through a material, and the absorbance A can be calculated as follows:

$$A = \lg (\frac{I_{v0}}{I_v})$$

(6)

where $I_{v0}$ is the luminous intensity passing through the empty Hele–Shaw cell, and $I_v$ is the luminous intensity after passing through the Hele–Shaw cell with solution.

As mentioned earlier, the driving circuitry of the table light is equipped with light brightness feedback to control the output light brightness. The light controller provides an adjustment function for 256 luminance levels. We measured the grayscale $G$ in the region of interest (ROI) under several different luminance levels $L_l$, concluded that the grayscale is linearly dependent on the luminance within a certain luminance range, and the relationship expressed as follows:

$$G =1.2L_l+45, L_l \in [ 28, 175 ].$$

(7)

The maximum luminance level within the linear correlation interval $L_l$ = 175 was selected as the incident light source luminance. In addition, luminance is a photometric measure of the luminous intensity per unit area of light traveling in a given direction, which is normally obtained by dividing the luminous intensity by the light source area.

Combining Equations (6) and (7), and for convenience in absorbance calculation, the ratio of luminous intensity can replace by the ratio of luminance level. The absorbance A can be translated as follows:

$$A=1\mathrm{g}\left(\frac{L_{l_0}}{L_l}\right)=-1\mathrm{g}\left(\frac{G - 45}{G_0 - 45}\right)$$

(8)

where $G_0$ is the gray value in the ROI of the empty Hele–Shaw cell image, and G is the gray value in the ROI of the Hele–Shaw cell with solution.

As the Beer–Lambert law [41] points out, absorbance is correlated with the concentrations of attenuating species as well as the thickness of the solution:

$$A = \alpha ·b·c$$

(9)

where α is the absorbance coefficient, b is the thickness of the solution, and c is the concentration of the absorbent solute.

The pH indicator consists of organic molecules that exhibit acid–base properties. The acidic form of the indicator (HIn) has a different color compared to the conjugate base (In⁻). A general equilibrium expression for the indicator is as follows:

$$\mathrm{HIn} \rightleftharpoons { \mathrm{H}}^{+}{+ \mathrm{In}}^{-}.$$

(10)

The equilibrium constant of Equation (10) is as follows:

$$K_{\mathrm{HIn}}{ = [\mathrm{H}}^{+}] \times \frac{{[\mathrm{In}}^{-}]}{\left[\mathrm{HIn}\right]}$$

(11)

where [ ] indicates the molar concentration of the species inside the square brackets.

The logarithm of Equation (11) can be simplified to the following:

$${\lg (K}_{\mathrm{HIn}}) = \lg \left(\frac{{[\mathrm{In}}^{-}]}{\left[\mathrm{HIn}\right]}\right)-\mathrm{pH}.$$

(12)

The hydrogen ion (H⁺) released due to the CO₂ dissolution changes the ratio of [HIn] to [${\mathrm{In}}^{-}$], which causes a color variation in the solution. The absorbance coefficient of ${\mathrm{In}}^{-}$ is ${\alpha }_{-}$, and the absorbance coefficient of HIn is α₀. The absorbance of the solution is as follows:

$$A = {\alpha }_{-}·b·\left[{\mathrm{In}}^{-}\right]+{\alpha }_0·b·\left[\mathrm{HIn}\right].$$

(13)

The pH indicator concentration c is equal to:

$$c = \left[{\mathrm{In}}^{-}\right] + \left[\mathrm{HIn}\right].$$

(14)

The absorbance of the indicator in a strong acid or strong alkali environment is as follows:

$$A_{\left(\mathrm{HIn}\right)}{ = \alpha }_0·b·\left[\mathrm{HIn}\right]{ = \alpha }_0·b·c$$

(15)

$$A_{{(\mathrm{In}}^{-})} = {\alpha }_{-}·b·\left[{\mathrm{In}}^{-}\right] ={\alpha }_{-}·b·c.$$

(16)

By using the extreme value of absorbance in Equation (13) and combining Equations (10), (11) and (14), the following equations are derived:

$$A = \frac{A_{{(\mathrm{In}}^{-})}}{c}{[\mathrm{In}}^{-}] + \frac{A_{\left(\mathrm{HIn}\right)}}{c}\left[\mathrm{HIn}\right]$$

(17)

$$\left[{\mathrm{In}}^{-}\right]\cdot \left(A_{\left({\mathrm{In}}^{-}\right)}-A\right) = \left[\mathrm{HIn}\right]\cdot (A-A_{\left(\mathrm{HIn}\right)})$$

(18)

$$\lg (\frac{K_{\mathrm{HIn}}}{{[\mathrm{H}}^{+}]}) =\lg (\frac{{[\mathrm{In}}^{-}]}{\left[\mathrm{HIn}\right]})=\lg \left(\frac{A-A_{\left(\mathrm{HIn}\right)}}{A_{\left({\mathrm{In}}^{-}\right)}-A}\right)$$

(19)

$$\lg \left(\frac{A-A_{\left(\mathrm{HIn}\right)}}{A_{\left({\mathrm{In}}^{-}\right)}-A}\right) = \mathrm{pH} -\mathrm{p}{K}_{HIn}$$

(20)

where $\lg \left(\frac{A-A_{\left(\mathrm{HIn}\right)}}{A_{\left({\mathrm{In}}^{-}\right)}-A}\right)$ is a factor related to absorbance and is proportional to the solution pH. A series of tests were conducted to calibrate the linear relationship between pH and $\lg \left(\frac{A-A_{\left(\mathrm{HIn}\right)}}{A_{\left({\mathrm{In}}^{-}\right)}-A}\right)$. Hydrochloric acid (HCl) or sodium hydroxide (NaOH) was added to the pH indicator solution to obtain five solutions with different pH values that were within the functional range of the indicator. The pH value was measured using a pH meter.

The extreme values of absorbance, A(HIn) and A(${\mathrm{In}}^{-}$), as shown in

Table 2. Extreme absorbance values and pK_HIn at different temperatures by experiments.

	25 °C	35 °C	45 °C
AHIn	0.039	0.044	0.051
AIn-	0.987	0.993	1.014
pKHIn	4.67	4.71	4.77

, were derived from the gray value of the images for the pH = 1 and pH = 10 indicator solutions by titration with HCl or NaOH.

As shown in Figure 2, there is good linearity between $\lg \left(\frac{A-A_{\left(\mathrm{HIn}\right)}}{A_{\left({\mathrm{In}}^{-}\right)}-A}\right)$ (the factor related to absorbance) and solution pH. This curve is the standard curve for the dissolved CO₂ concentration measurement experiments, which is used to calculate [H⁺] from the absorbance.

Based on the requirement for electrical neutrality in the solution, the following charge balance equation is obtained:

$$\left[{\mathrm{HCO}}_3^{-}\right] + \left[{\mathrm{In}}^{-}\right]{ = [\mathrm{H}}^{+}]+ \left[{\mathrm{Na}}^{+}\right]$$

(21)

where the $\left[{\mathrm{Na}}^{+}\right]$ is equal to the concentration of sodium hydroxide in the pH indicator solution.

Using the equilibrium constant K_HIn (Equation (11)) and the indicator concentration c = 2.5 × 10⁻⁴ mol/L from Equation (18), [In^-] can be written as follows:

$$\left[{\mathrm{In}}^{-}\right] = \frac{c}{{10}^{-\mathrm{pH}}·{ K}_{\mathrm{HIn}}^{-1} + 1}.$$

(22)

By combining the standard curve with Equations (17)–(19), ${[\mathrm{HCO}}_3^{-}]$ can be calculated.

The equilibrium expression involving CO₂ (aq) (Equation (5)) at different temperatures [43] can be combined with the equilibrium constant for the reaction and written as follows:

$$\left[{\mathrm{CO}}_2\left(\mathrm{aq}\right)\right] = \frac{\left[{\mathrm{H}}^{+}\right] ·{ [\mathrm{HCO}}_3^{-}]}{K_{{\mathrm{CO}}_2\left(\mathrm{aq}\right)}}.$$

(23)

Finally, the concentration of dissolved CO₂ consists of two parts: [CO₂(aq)] and [${\mathrm{HCO}}_3^{-}$]:

$${[\mathrm{CO}}_2{\left(\mathrm{dissolved}\right)] = [\mathrm{CO}}_2{\left(\mathrm{aq}\right)] + [\mathrm{HCO}}_3^{-}].$$

(24)

3.3. Numerical Simulation

The numerical model used in our simulation was built similar to the experimental configuration considering the CO₂-saturated water interface in the Hele–Shaw cell. The top boundary is at atmospheric pressure conditions. Bottom and sides boundaries are considered as no-flow boundaries. Simulation was performed using STOMP (Subsurface Transport Over Multiple Phases), and the operational mode is STOMP-CO2.

The governing equations of flow and concentration field for such a system are as follows:

$$\frac{\partial }{\partial t}\left[{{\displaystyle \sum }}^{\hspace{0.17em}}\left(\rho \varphi S{\omega }^i\right)\right]=-{{\displaystyle \sum }}^{\hspace{0.17em}}\nabla \left(\rho {\omega }^{i}V\right)-{{\displaystyle \sum }}^{\hspace{0.17em}}\nabla \left(F^i\right)+{{\displaystyle \sum }}^{\hspace{0.17em}}\left(\omega m_d\right)$$

(25)

where S is the saturation, ω is the mass fraction, i represents the species (H₂O or CO₂), and m_d is the mass rate density.

Darcy’s law is used to compute the advective fluxes V:

$$\mathit{V}=-\frac{k_{r}k}{\mu }\left(\nabla P+\rho gh\right).$$

(26)

The diffusive flux F is calculated considering molecular diffusion and neglect dispersion:

$${\mathit{F}}^i=- \rho \varphi S\frac{M^i}{M}\left(\tau D^i\right)\nabla {\mathsf{\chi }}^i$$

(27)

where M is molecular weight, $\tau$ is tortuosity, and $\chi$ is mole fraction.

The saturation–capillary pressure function is that of van Genuchten [44]:

$$\overline{S_{el}}={\left[1+{\left(\alpha h\right)}^n\right]}^{-m}; m=1-\frac{1}{n}; \overline{S_{el}}=\frac{S_l-S_r}{1-S_r}.$$

(28)

The relative permeability function of van Genuchten ($\alpha$ = 500; n = 5) is associated with the Mualem porosity distribution model [45,46]:

$$k_r=\sqrt{\overline{S_{el}}}\left[{\left(1-{\left(1-{\left(\overline{S_{el}}\right)}^{1/m}\right)}^m\right)}^2\right].$$

(29)

A detailed description of all the state equations used in STOMP-CO2 can be found in White et al. [47]. The mesh sensitivity analysis was conducted, and we concluded that the uniform mesh of 100 × 100 (grid cell length and height equal to 1 mm) was adequate when the CO₂ mass transfer remained stable. The simulation was perturbed with non-regular sinusoidal perturbation introduced on the initial concentration profile just below the interface.

4. Results and Discussion

The main aim of this work was to establish a new spectrophotometric method for the quantitative visualization of CO₂ convective dissolution in brine, investigate the effects of salinity and temperature on convective instability and the mass of CO₂ dissolved, and finally use the experimental data for simulation validation.

4.1. Development Regimes and Morphology of Instabilities

Figure 3 shows successive images of the CO₂ dissolution in deionized water at 25 °C. The selected ROI has a dimension of 75 mm × 75 mm. When CO₂ was injected into the Hele–Shaw cell with brine, CO₂ started to dissolve into the brine, and the CO₂-dissolved brine was denser than the native brine. The CO₂-dissolved brine, which is shown as a thin and bright line layer on the images (Figure 3a), appeared almost instantaneously behind the interface once the CO₂ was injected into the Hele–Shaw cell, and this indicated the beginning of the reaction. At the initial stage, the only active mechanism is diffusion, and the denser brine layer below the interface is relatively stable; this stage refers to the diffusion dissolution regime. Next, the CO₂ mass influx produced an increase in density, the density difference caused instability, small fingers formed, and eventually, convection was triggered.

At early times, convective fingers grow independently (Figure 3b). As time passes, there are interactions between adjacent fingers. The small fingers first merged at the interface, and then the tips of the fingers tended to coalesce into longer and larger fingers with a few seconds of elapsed time. The wave numbers were determined by counting the number of fingers at the top (interface) and the bottom (tip of fingers) using the length scale. It is clearly seen from 600 to 11,000 s (Figure 3c–k) that as the process proceeds, convective fingers are fully developed, and wave numbers decreased both at the top and bottom. The CO₂-dissolved brine moved downwards, and the native brine moved upwards to the interface. Meanwhile, new fingers emerged from the interface between existing fingers and follow the same behavior. Several dynamic regimes of CO₂ convective dissolution have been proposed in the literature [37,48]; this period in our experiments is known as the convective regime. Eventually, the downward fingers touched the bottom boundary and became blurred, and this resulted in the decay of the dissolution as the system entered the shutdown regime (Figure 3k–l).

When the convection regime is active, the dissolution rate of CO₂ into the brine is greater. As can be seen from Figure 3, within 100 s after introducing CO₂ into the top of the brine, only a thin-layer area below the interface became bright, and diffusion dominates the mass transport. However, from time t = 600 s, when convection dominates the mass transport, a larger area of brine became bright, which shows the direct effects of convection on the enhanced CO₂ dissolution. The development regimes variation of density-driven convective are comparable to the reported results [21]. We will perform a quantitative analysis of CO₂ dissolution in the later section.

4.2. Effects of Brine Salinity and Temperature

To investigate the effects of salinity on convective instability of CO₂ dissolution, we analyzed the convection starting time and the length of fingers. The convection starting time is the time at which the convective fingers can be observed. With the careful observation of our experimental images, the approximate convection starting time (based on the first appearance of fingers at the interface) is recorded. Under the experimental conditions in this study, the convection starting time occurred within 300 s. Analytical and numerical calculations indicated that the convection starting time shows a relation with $\approx 1/R{a}^2$ [49]. We plotted the average value of the convection starting time change with Rayleigh number, as shown in Figure 4. The same behavior is observed in the variation of the convection starting time in the different salinities by the experimental measurement and linear stability analysis calculation. The higher the salinity, the smaller the Ra and the longer the convection starting time. In conclusion, increasing salinity changes the convection starting time in an undesirable way.

Figure 5 shows the CO₂ convective dissolution images at 3000 s for different salinities. The length of the fingers was measured from the interface to the tip of the fingers. These four experiments exhibit different instability behavior in different fingering pattern. Lower salinity results in a faster development of instabilities such that the growth rate of fingers is significantly stronger. The temporal evolution of the fingers length in different salinity conditions is shown in Figure 6, where the observed growth is not linear. The change in fingers length at the initial times is fast compared with final times. After the short diffusion stage, a sharp increase indicates the evolution of a convective process.

When CO₂ fingers occurred at the highest salinity (Exp. No. 4), compared with other conditions, it was the shortest growing fingers. The variation between the three curves in Figure 6 is due to the differences in density caused by salinity. The CO₂ solubility in brine decreased with increasing salinity, and the density difference between brine with and without dissolved CO₂ is decreased [50]. The viscosity of aqueous solutions also increased with increasing salinity because it is based on the density variation. Dissolved ionic compounds increase the density and therefore its viscosity [21]. The density difference is the major factor that determines how fast the convective fingers move [51]. In short, salinity has a negative effect on CO₂ convective dissolution. With increasing salinity, the length of fingers decreases, and the fingers growth rate decreases.

Our experiments were also performed under different temperature conditions, varying the temperature from 25 to 45 °C. The increase in temperature led to the increased Rayleigh number, and the appearance time of the largest wave number is delayed, as shown in Figure 7. The CO₂-affected area of the convective pattern was analyzed for different temperature conditions. This area corresponds to the ratio of the fingers area to the total aqueous phase area and could be interpreted as an approximate measure of CO₂ dissolution. The larger the affected area, the larger the amount of dissolved CO₂. To distinguish the CO₂-affected area from the area without dissolved CO₂, the images were converted into binary images through threshold segmentation, and the affected area can be clearly shown on the background.

To make a qualitative comparison, the ImageJ area fraction measurement function was applied to obtain the percentage of the CO₂-affected area. The affected areas over different temperature experiments are plotted as a function of time in Figure 8. With CO₂ dissolution, the affected area rapidly increased and then gradually stabilized. It should be noted that due to the temperature difference between these experiments, the CO₂-affected area in deionized water at 25 °C was larger than the others, and the final values of the CO₂-affected area decreased with increasing temperature.

Based on thermodynamic theory, the apparent molar volume, that is one parameter of the CO₂ aqueous solution density, is a function of temperature [52,53]. The increase in the experimental temperature leads to a decrease in the CO₂ aqueous solution density, which further causes the density difference to decrease. The density difference in the experimental fluids is detailed in Table 1. The density difference is the driving force for convection, and the density difference caused by dissolved CO₂ increases the CO₂-affected area. In addition, as CO₂ reacts with water, its solubility depends on the temperature condition, and CO₂ solubility decreases with temperature under isobaric conditions [54]. The lower the solubility, the less amount of CO₂ that could react with water. From the curves, it can be concluded that with increasing temperature, the CO₂-affected area decreases, and the amount of dissolved CO₂ decreases.

4.3. Determination of the Mass of Dissolved CO₂

To determine and compare the concentration of dissolved CO₂ under all experimental conditions, we selected the first 15,000 s (250 min) of the experiment as the characteristic time. Before running the actual experiments, reliability tests were performed using HCl and NaOH solutions of known pH values to verify the accuracy of the method. For three different pH values, the experimental error was found to be less than 5%. In addition, we also measured the solution absorbance with different indicator concentrations and calculated the mass absorption coefficient of the pH indicator at the maximum wavelength. The experimental value has a small error (3.8%) to the property parameter. It has also been suggested that the experimental system is reliable for measuring absorbance.

As described in Section 3.2, the cumulative dissolved CO₂ concentration was determined by the gray value of the images combined with spectrophotometric analysis. Figure 9 shows the change in the cumulative dissolved CO₂ concentration for different experimental temperatures and salinities. The dissolved CO₂ concentration in Exp. No. 1 is greater than all other curves, and this is due to the decrease in salinity and temperature.

The solubility of Experiment (Exp.) No. 1 at 15,000 s is 1.34 kg/m³, and the errors between our experimental result and those of the solubility model of Duan et al. [39] is 6.3%. However, under other experimental conditions, the solubility of our experiments at 15,000 s deviates from the model because the calculated data for convection are not fully developed (see Figure 10); as the bright area of the images increases, the amount of dissolved CO₂ becomes larger.

After calculating the cumulative concentration of dissolved CO₂, the cumulative mass of dissolved CO₂ is calculated. Fick’s second law and the following equation, presented by Amir Taheri et al. [55], are used to calculate the dissolved CO₂ mass of diffusion dissolution:

$$M\left(t\right)=2C_0\sqrt{\frac{Dt}{\pi }}.$$

(30)

In the experimental setup suggested here, the concentration is maintained at a constant C₀. Here, the CO₂ equilibrium concentration at the interface was determined by the solubility model [39].

Reducing the diffusion dissolution of CO₂ from cumulative dissolved CO₂ results in the convective dissolution of CO₂. Figure 11 displays the dissolved mass of CO₂ of Exp. No. 1 for different dissolution mechanisms. By comparing the mass of diffusion dissolved CO₂ to convection dissolved CO₂, the crucial role of convection on the enhancement of CO₂ mass transfer can be seen. The curve also indicates that the convective dissolution effect increases within 15,000 s. By differentiating from the mass of cumulative dissolved CO₂, the dissolved CO₂ flux for each experiment is reached. This investigation is presented in Figure 12.

The CO₂ dissolution flux curves have different maximum values at the initial convection regime, and they decrease while increases in temperature and salinity cause a decrease in the dissolution flux. Earlier studies on convective systems have focused on the Sherwood number and Rayleigh number. Under our experimental conditions and by increasing salinity and temperature, variables such as density differences and solubility would decrease, the viscosity and diffusion coefficient would increase; thus, their total effects reduce the Rayleigh number.

To find a relationship between Ra and Sh, the maximum dissolution flux was used. By plotting Sh against Ra, a relationship is obtained as shown in Figure 13. Here, Sh scales with Ra as $Sh=0.203R{a}^{0.832}$. This agrees with existing numerical [28,56] and experimental [21,31] results with Rayleigh number exponents approximately equal to 4/5. These slight variations in the value of the exponent are attributed to the error in calculating the CO₂ dissolution flux, the differences in fluid systems and the experimental setups. Figure 13 displays that Sh increases with increasing Ra, which means that the mass transfer flux of CO₂ dissolution increases with salinity and as the temperature decreases. If the physical parameters of the candidate field are known, the dimensionless CO₂ dissolution flux can be estimated, thereby helping to examine the strength of geological CO₂ storage.

4.4. Comparison between Experimental and Simulation Results

The numerical simulation uses the experimental results for comparison and validation of the patterns of convective fingers and the mass of dissolved CO₂. Figure 14 shows the morphology of the fingers in experiment and simulation for the CO₂ dissolution in deionized water at 25 °C. The simulation results and experimental results accord with very well considering the fingers pattern at the same times, and the fingers become longer and thicker over time. Compared with the experimental results, the simulation show that the rate of fingers growth is slightly faster. It is because the boundary condition of constant saturation in dissolved CO₂ at the top of the Hele–Shaw cell shorten the time for CO₂ to diffuse into the aqueous phase.

Figure 15 simulates Experiments 2–6 of different salinity and temperature conditions at 3000 s. The variation trend of fingers length and influence area in the simulation is similar to the experimental results. The increase in brine salinity and temperature adversely affects the growth and development of the fingers. Moreover, the average position of the fingertip in the simulation is lower in comparison to the experiments (see Figure 5).

Figure 16 compares the total mass of dissolved CO₂ at 3000 s for different conditions. The order of magnitudes for the dissolved CO₂ mass in the simulation is the same as that in the experiments. However, there are some differences between the simulation and experimental results for those salinity and temperature variations cases. The numerical simulations overestimated the dissolved CO₂ mass affected by salinity and temperature. Comparison of the experimental results with the numerical results shows that the dissolved CO₂ mass decreases with the increase of salinity and temperature. Moreover, the regularity of the dissolved CO₂ mass in the simulation (red lines) is more than the experiments (black lines). The assumptions used for the simulation interface seemed to be the primary cause for the discrepancy. The main similarity between experimental and simulation results is that CO₂ convective dissolution is more effective for CO₂ storage reservoir with lower salinity and temperature.

5. Conclusions

This study developed a quantitative method for analyzing the density-driven convective instability and measuring the mass transfer of CO₂ dissolution. Several experiments and simulation were performed to investigate the effects of brine salinity and temperature on CO₂ dissolution relevant to CO₂ storage in saline aquifers. Based on the analysis of the results, the following conclusions can be drawn:

• There are three different regimes of instability development: diffusion dissolution regime, convection diffusion regime, and shutdown regime. The course and development of the convective fingers go through generation, propagation, coalescence, and re-initiation.
• Brine salinity has a negative effect on the CO₂ convective dissolution. Thus, with increases in salinity, the convection starting time delay, the length of fingers, and fingers growth rate decrease.
• Although increasing the temperature causes an increase in the diffusion coefficient, the solubility of CO₂ in water decreases the density difference. An increase in temperature leads to a reduction in the initial wave number and CO₂ affected area, which can adversely affect CO₂ dissolution.
• The spectrophotometric method was successfully used to quantitatively measure the dissolved CO₂ concentration. The mass of dissolved CO₂ due to convective dissolution is much greater than that for diffusion dissolution. The dimensionless flux, where the Sherwood number shows a power-law relationship with the Rayleigh number, indicates that the Sh increases with increasing Ra, and this refers to the lower salinity and temperature increasing the mass transfer flux of CO₂ dissolution.
• Salinity and temperature effects on convective instability of numerical simulation are in qualitative agreement with the experimental result. The shape and position of the convective finger in the simulation are influenced by the assumption of boundary conditions. Numerical simulation of CO₂ convective dissolution still requires further research.

A quantitative and visualization study of the CO₂ dissolution process is important for the further development of technology and simulation models. A determination of the dissolution flux can help to estimate the amount of CO₂ that will remain in brine during the CO₂ storage in saline aquifers. Our findings improve the understanding of the CO₂ dissolution mechanism and help to assess the CO₂ solubility trapping rate.