Convective dissolution of carbon dioxide in deep saline aquifers: Insights from engineering a high-pressure porous Hele-Shaw cell

TL;DR

This study experimentally investigates CO2 dissolution-driven convection in high-pressure porous media using a novel Hele-Shaw cell, revealing dynamic regimes, feedback effects, and scaling relations relevant for geologic carbon storage.

Contribution

First experimental analysis of high-pressure CO2 convection in porous media with a new Hele-Shaw cell, providing insights into regimes, feedback mechanisms, and scaling laws.

Findings

Dissolution flux decreases continuously in the intermediate regime.

NaCl+CaCl2 mixtures extend the intermediate period but reduce flux.

Introduced a new scaling factor for flux and regime transition times.

Abstract

We present the first experiments of dissolution-driven convection of carbon dioxide (CO$_2$) in a confined brine-saturated porous medium at high pressures. We designed a novel Hele-Shaw cell that allows for both visual and quantitative analyses, and address the effects of free-phase CO$_2$ and brine composition on convective dissolution. The visual examination of the gas volume combined with the measurement of pressure, which both evolve with dissolution, enable us to yield insights into the dynamics of convection in conditions that more closely reflect the geologic conditions. We find and analyze different dissolution events, including diffusive, early and late convection, and shut-down regimes. Our experiments reveal that in intermediate regime, a so-called `quasi-steady' state actually never happens. Dissolution flux continuously decreases in this regime, which is due to a negative feedback loop: the rapid reduction of pressure following convective dissolution, in turn, decreases the solubility of CO$_2$ at the gas-brine interface and thus the instability strength. We introduce a new scaling factor that not only compensates the flux reduction but also the nonlinearities that arise from different salt types. We present robust scaling relations for the compensated flux and for the transition times between consecutive regimes in systems with NaCl (Ra $\sim$ 3271--4841) and NaCl+CaCl$_2$ mixtures (Ra $\sim$ 2919--4283). We also find that NaCl+CaCl$_2$ mixtures enjoy a longer intermediate period before the shut-down of dissolution, but with a lower dissolution flux, as compared to NaCl brines...

Tables

Table 1. Table 1: A summary of scaling relations for estimating dissolution flux. ∗ ∗ \ast Methanol and ethylene-glycol     ∗ ∗ \ast ∗ ∗ \ast Propylene glycol     ∗ ∗ \ast ∗ ∗ \ast ∗ ∗ \ast a, a ′ are constant values

Reference	System	Regime	Analysis method	Ra range	Equation
\cellcolor[HTML]FFFFFF(Hesse, 2008)	\cellcolor[HTML]FFFFFF-	\cellcolor[HTML]FFFFFFQuasi-steady	\cellcolor[HTML]FFFFFFNumerical	\cellcolor[HTML]FFFFFF-	\cellcolor[HTML]FFFFFF $Sh=0.017Ra$
\cellcolor[HTML]FFFFFF(Neufeld et al., 2010)	\cellcolor[HTML]FFFFFFMEG$\ast$ -water	\cellcolor[HTML]FFFFFFQuasi-steady	\cellcolor[HTML]FFFFFFExperimental Numerical	\cellcolor[HTML]FFFFFF $2\times {10}^3-$ $6\times {10}^5$	\cellcolor[HTML]FFFFFF $Sh=0.12R{a}^{0.84}$
\cellcolor[HTML]FFFFFF(Pau et al., 2010)	\cellcolor[HTML]FFFFFF-	\cellcolor[HTML]FFFFFFQuasi-steady	\cellcolor[HTML]FFFFFFNumerical	\cellcolor[HTML]FFFFFF-	\cellcolor[HTML]FFFFFF $F=0.017\frac{c_{0}k\Delta \rho g}{\mu }$
\cellcolor[HTML]FFFFFF(Backhaus et al., 2011)	\cellcolor[HTML]FFFFFFPPG$\ast$ $\ast$ -water	\cellcolor[HTML]FFFFFFQuasi-steady	\cellcolor[HTML]FFFFFFExperimental	\cellcolor[HTML]FFFFFF $6\times {10}^3-$ $9\times {10}^4$	\cellcolor[HTML]FFFFFF $Sh=0.045R{a}^{0.76}$
\cellcolor[HTML]FFFFFF(Elenius and Johannsen, 2012)	\cellcolor[HTML]FFFFFF-	\cellcolor[HTML]FFFFFFQuasi-steady	\cellcolor[HTML]FFFFFFNumerical	\cellcolor[HTML]FFFFFF-	\cellcolor[HTML]FFFFFF $F=0.02\frac{c_{0}k\Delta \rho g}{\mu }$
\cellcolor[HTML]FFFFFF(Hidalgo et al., 2012)	\cellcolor[HTML]FFFFFF-	\cellcolor[HTML]FFFFFFQuasi-steady	\cellcolor[HTML]FFFFFFNumerical	\cellcolor[HTML]FFFFFF $5\times {10}^3-$ $3\times {10}^4$	\cellcolor[HTML]FFFFFF $F=a\frac{c_{0}k\Delta \rho g}{\mu }$
\cellcolor[HTML]FFFFFF(Farajzadeh et al., 2013)	\cellcolor[HTML]FFFFFF-	\cellcolor[HTML]FFFFFFQuasi-steady	\cellcolor[HTML]FFFFFFNumerical	\cellcolor[HTML]FFFFFF $1\times {10}^3-$ $8\times {10}^3$	\cellcolor[HTML]FFFFFF $Sh=0.0794R{a}^{0.832}$
\cellcolor[HTML]FFFFFF(Tsai et al., 2013)	\cellcolor[HTML]FFFFFFPPG$\ast$ $\ast$ -water	\cellcolor[HTML]FFFFFFQuasi-steady	\cellcolor[HTML]FFFFFFExperimental	\cellcolor[HTML]FFFFFF $5\times {10}^3-$ $1\times {10}^5$	\cellcolor[HTML]FFFFFF $Sh=0.037R{a}^{0.84}$
\cellcolor[HTML]FFFFFF(Elenius et al., 2014)	\cellcolor[HTML]FFFFFF-	\cellcolor[HTML]FFFFFFQuasi-steady	\cellcolor[HTML]FFFFFFNumerical	\cellcolor[HTML]FFFFFF-	\cellcolor[HTML]FFFFFF $F=0.021\frac{c_{0}k\Delta \rho g}{\mu }$
\cellcolor[HTML]FFFFFF(Green and Ennis-King, 2014)	\cellcolor[HTML]FFFFFF-	\cellcolor[HTML]FFFFFFQuasi-steady	\cellcolor[HTML]FFFFFFNumerical	\cellcolor[HTML]FFFFFF-	\cellcolor[HTML]FFFFFF $F=0.017\sqrt{k_v{k}_h}\frac{c_0\Delta \rho g}{\mu }$
\cellcolor[HTML]FFFFFF(Mojtaba et al., 2014)	\cellcolor[HTML]FFFFFFCO2-brine (NaCl)	\cellcolor[HTML]FFFFFFQuasi-steady	\cellcolor[HTML]FFFFFFExperimental	\cellcolor[HTML]FFFFFF182-20860	\cellcolor[HTML]FFFFFF $Sh=0.0228R{a}^{0.7897}$
\cellcolor[HTML]FFFFFF(Slim, 2014)	\cellcolor[HTML]FFFFFF-	\cellcolor[HTML]FFFFFFQuasi-steady	\cellcolor[HTML]FFFFFFNumerical	\cellcolor[HTML]FFFFFF $2\times {10}^3-$ $5\times {10}^5$	\cellcolor[HTML]FFFFFF $F=0.017\frac{c_{0}k\Delta \rho g}{\mu }$
\cellcolor[HTML]FFFFFF	\cellcolor[HTML]FFFFFF	\cellcolor[HTML]FFFFFFShut-down	\cellcolor[HTML]FFFFFF	\cellcolor[HTML]FFFFFF	\cellcolor[HTML]FFFFFF $Sh=\frac{16.8}{{\left[0.73\left(\frac{t}{Ra}-16\right)+31.5\right]}^2}$
\cellcolor[HTML]FFFFFF(Martinez and Hesse, 2016)	\cellcolor[HTML]FFFFFF-	\cellcolor[HTML]FFFFFFQuasi-steady	\cellcolor[HTML]FFFFFFNumerical	\cellcolor[HTML]FFFFFF-	\cellcolor[HTML]FFFFFF $F=(0.018,0.019)\frac{c_{0}k\Delta \rho g}{\mu }$
\cellcolor[HTML]FFFFFF(Wang et al., 2016)	\cellcolor[HTML]FFFFFFNaCl solution- MEG	\cellcolor[HTML]FFFFFFQuasi-steady	\cellcolor[HTML]FFFFFFExperimental	\cellcolor[HTML]FFFFFF2600-16036	\cellcolor[HTML]FFFFFF $Sh=0.13R{a}^{0.93}$
\cellcolor[HTML]FFFFFF(De Paoli et al., 2017)	\cellcolor[HTML]FFFFFF- \cellcolor[HTML]FFFFFF	\cellcolor[HTML]FFFFFFQuasi-steady	\cellcolor[HTML]FFFFFFNumerical \cellcolor[HTML]FFFFFF	\cellcolor[HTML]FFFFFF $1\times {10}^3-$ \cellcolor[HTML]FFFFFF $2\times {10}^4$	\cellcolor[HTML]FFFFFF $F=0.017\sqrt{k_v{k}_h}\frac{c_0\Delta \rho g}{\mu }$
\cellcolor[HTML]FFFFFF	\cellcolor[HTML]FFFFFF	\cellcolor[HTML]FFFFFFShut-down	\cellcolor[HTML]FFFFFF	\cellcolor[HTML]FFFFFF	\cellcolor[HTML]FFFFFF $Sh=\frac{4a{\left(\frac{k_v}{k_h}\right)}^{a^{\prime }}}{{\left(1+4a{\left(\frac{k_v}{k_h}\right)}^{a^{\prime }}t/Ra\right)}^2}$
\cellcolor[HTML]FFFFFF(Taheri et al., 2017)	\cellcolor[HTML]FFFFFFCO2-water	\cellcolor[HTML]FFFFFFQuasi-steady	\cellcolor[HTML]FFFFFFExperimental	\cellcolor[HTML]FFFFFF709- 9627	\cellcolor[HTML]FFFFFF $F=\frac{0.021\Delta \rho g.cos\left(\theta \right)k{c}_0}{\mu }$
\cellcolor[HTML]FFFFFF(Newell et al., 2018) \cellcolor[HTML]FFFFFF	\cellcolor[HTML]FFFFFFCO2-water \cellcolor[HTML]FFFFFF	\cellcolor[HTML]FFFFFFFlux-growth	\cellcolor[HTML]FFFFFFExperimental \cellcolor[HTML]FFFFFFNumerical	\cellcolor[HTML]FFFFFF2093-16256 \cellcolor[HTML]FFFFFF	\cellcolor[HTML]FFFFFF $F=a^{\prime }{t}^{0.5}$
\cellcolor[HTML]FFFFFF	\cellcolor[HTML]FFFFFF	\cellcolor[HTML]FFFFFFQuasi-steady	\cellcolor[HTML]FFFFFF	\cellcolor[HTML]FFFFFF	\cellcolor[HTML]FFFFFF $F=a\frac{c_{0}k\Delta \rho g}{\mu \phi }$
\cellcolor[HTML]FFFFFF	\cellcolor[HTML]FFFFFF	\cellcolor[HTML]FFFFFFShut-down	\cellcolor[HTML]FFFFFF	\cellcolor[HTML]FFFFFF	\cellcolor[HTML]FFFFFF $F=a^{{}^{\prime }}{t}^{-1.75}$
\cellcolor[HTML]FFFFFF(Wen et al., 2018)	\cellcolor[HTML]FFFFFF-	\cellcolor[HTML]FFFFFFQuasi-steady	\cellcolor[HTML]FFFFFFNumerical	\cellcolor[HTML]FFFFFF $1\times {10}^4-$ $5\times {10}^4$	\cellcolor[HTML]FFFFFF $Sh=\frac{0.0168R{a}_0}{{\left(0.0168a{t}_a+1\right)}^2}$
\cellcolor[HTML]FFFFFF	\cellcolor[HTML]FFFFFF	\cellcolor[HTML]FFFFFFShut-down	\cellcolor[HTML]FFFFFF	\cellcolor[HTML]FFFFFF	\cellcolor[HTML]FFFFFF $Sh=\frac{0.0317R{a}_0}{{\left(0.0317\left(1+a\right)t_a+0.861\right)}^2}$

Table 2. Table 2: The details of experimental tests with their order number, brine composition, permeability, Rayleigh value, initial pressure, and diffusion coefficient. The details of correlations that are used to calculate Ra are provided in Appendix A.

Test Number	brine composition (mole based)	Permeability (D)	Rayleigh number	Initial Pressure (psi)	Diffusion $\times$${10}^9\left(\frac{m^2}{s}\right)$
1	2 NaCl	550	4444	523.1	3.8
2	2 NaCl	400	3272	511.2	3.7
3	1 NaCl	550	4841	535.3	5
4	1 NaCl	400	3514	510.6	4.9
5	1.6 NaCl+ 0.2 CaCl2	550	3893	512.7	4.2
6	1.6 NaCl + 0.2 CaCl2	400	2919	514.7	4.1
7	0.8 NaCl + 0.1 CaCl2	550	4283	502.6	5.3
8	0.8 NaCl + 0.1 CaCl2	400	3242	505.5	5.2

Table 3. Table 3: Fitted and estimated factors for the modified dissolution flux equation (after (Wen  et al. , 2018 ) ) and transition times between different dissolution regimes

Case	$\frac{F_c}{C_s^2}$	$m_{fit.}^{F.G}$ $\times {10}^4$	$b_{fit.}^{F.G}$ $\times {10}^6$	$m_{est.}^{F.G}$ $\times {10}^4$	$b_{est.}^{F.G}$ $\times {10}^6$	$t_{D,fit.}^{Q.S}$ $\times {10}^4$	$t_{D,est.}^{Q.S}$ $\times {10}^4$	$a_{fit.}^{Q.S}$ $\times {10}^4$	$a_{est.}^{Q.S}$ $\times {10}^4$	$t_{D,fit.}^{Sh}$ $\times {10}^4$	$t_{D,est.}^{Sh}$ $\times {10}^4$	$a_{fit.}^{S.D}$ $\times {10}^9$
1	7	77	-124	78	-134	28.73	28.58	2.90	2.95	54.13	53.03	8.3
\cellcolor[HTML]FFFFFF2	\cellcolor[HTML]FFFFFF5	\cellcolor[HTML]FFFFFF45	\cellcolor[HTML]FFFFFF-52	\cellcolor[HTML]FFFFFF43	\cellcolor[HTML]FFFFFF-50	\cellcolor[HTML]FFFFFF38.47	\cellcolor[HTML]FFFFFF37.64	\cellcolor[HTML]FFFFFF2.20	\cellcolor[HTML]FFFFFF2.09	\cellcolor[HTML]FFFFFF71.00	\cellcolor[HTML]FFFFFF70.47	\cellcolor[HTML]FFFFFF11
3	4.1	96	-173	90	-162	26.47	26.46	3.30	3.24	48.10	48.99	6.5
\cellcolor[HTML]FFFFFF4	\cellcolor[HTML]FFFFFF3.4	\cellcolor[HTML]FFFFFF55	\cellcolor[HTML]FFFFFF-65	\cellcolor[HTML]FFFFFF50	\cellcolor[HTML]FFFFFF-68	\cellcolor[HTML]FFFFFF34.27	\cellcolor[HTML]FFFFFF35.30	\cellcolor[HTML]FFFFFF2.30	\cellcolor[HTML]FFFFFF2.26	\cellcolor[HTML]FFFFFF64.47	\cellcolor[HTML]FFFFFF65.95	\cellcolor[HTML]FFFFFF5.5
5	4.5	61	-92	62	-95	31.76	32.14	2.40	2.54	85.76	85.88	18.0
\cellcolor[HTML]FFFFFF6	\cellcolor[HTML]FFFFFF3.8	\cellcolor[HTML]FFFFFF31	\cellcolor[HTML]FFFFFF-31	\cellcolor[HTML]FFFFFF32	\cellcolor[HTML]FFFFFF-25	\cellcolor[HTML]FFFFFF42.30	\cellcolor[HTML]FFFFFF42.99	\cellcolor[HTML]FFFFFF1.70	\cellcolor[HTML]FFFFFF1.84	\cellcolor[HTML]FFFFFF104.9	\cellcolor[HTML]FFFFFF105.18	\cellcolor[HTML]FFFFFF18.2
7	4.6	64	-121	73	-123	29.10	29.19	2.80	2.83	80.31	80.30	17
\cellcolor[HTML]FFFFFF8	\cellcolor[HTML]FFFFFF3.2	\cellcolor[HTML]FFFFFF42	\cellcolor[HTML]FFFFFF-46	\cellcolor[HTML]FFFFFF42	\cellcolor[HTML]FFFFFF-48	\cellcolor[HTML]FFFFFF39.32	\cellcolor[HTML]FFFFFF38.67	\cellcolor[HTML]FFFFFF2.20	\cellcolor[HTML]FFFFFF2.07	\cellcolor[HTML]FFFFFF98.14	\cellcolor[HTML]FFFFFF97.69	\cellcolor[HTML]FFFFFF15.6

Table 4. Table 4: Parameters of c 1 -c 15 , required in CO 2 solubility equations

Constant	$\frac{{\mu }_{{\mathrm{CO}}_2}^{1\left(0\right)}}{RT}$	${\lambda }_{{\mathrm{CO}}_2-\mathrm{Na}}$	${\zeta }_{{\mathrm{CO}}_2-\mathrm{Na}-\mathrm{Cl}}{m}_{\mathrm{Cl}}$	${\varphi }_{{\mathrm{CO}}_2}$
$c_1$	28.9447706	-0.411370585	3.36389723E-4	1.0
$c_2$	-0.0354581768	6.07632013E-4	-1.98298980E-5	4.7586835E-3
$c_3$	-4770.67077	97.5347708	-	-3.3569963E-6
$c_4$	1.02782768E-5	-	-	0.0
$c_5$	33.8126098	-	-	-1.3179396
$c_6$	9.04037140E-3	-	-	-3.8389101E-6
$c_7$	-1.14934031E-3	-	-	0.0
$c_8$	-0.307405726	-0.0237622469	2.12220830E-3	2.2815104E-3
$c_9$	-0.0907301486	0.0170656236	-5.24873303E-3	0.0
$c_{10}$	9.32713393E-4	-	-	0.0
$c_{11}$	-	1.41335834E-5	-	0.0

Table 5. Table 5: Required parameters for the viscosity models

Constant	Mao Duan’s model (Mao and Duan, 2009)	Zhang et al. model (Zhang et al., 1997)
$c_1$	-0.21319213	0.0061
$c_2$	0.13651589E-2	-
$c_3$	-0.12191756E-5	0.01040
$c_4$	0.69161945E-1	0.000756
$c_5$	-0.27292263E-3	-
$c_6$	0.20852448E-6	0.0157
$c_7$	-0.25988855E-2	0.271
$c_8$	0.77989227E-5	0.04712
$c_9$	-	0.00941
$c_{10}$	-	0.00003

Table 6. Table 6: Parameters for the calculation of NaCl + CaCl 2 brine density

\cellcolor[HTML]FFFFFFCon.	\cellcolor[HTML]FFFFFF${\alpha }_{10}$	\cellcolor[HTML]FFFFFF${\alpha }_{11}$	\cellcolor[HTML]FFFFFF${\alpha }_{12}$	\cellcolor[HTML]FFFFFF${\alpha }_{13}$	\cellcolor[HTML]FFFFFF${\alpha }_{14}$	\cellcolor[HTML]FFFFFF${\alpha }_{20}$	\cellcolor[HTML]FFFFFF${\alpha }_{21}$
NaCl	2863.158	-46844.356	120760.118	-116867.722	40285.426	-2000.028	34013.518
CaCl2	2546.760	-39884.946	102056.957	-98403.334	33976.048	-1362.157	22785.572
\cellcolor[HTML]FFFFFFCon.	\cellcolor[HTML]FFFFFF${\alpha }_{22}$	\cellcolor[HTML]FFFFFF${\alpha }_{23}$	\cellcolor[HTML]FFFFFF${\alpha }_{24}$	\cellcolor[HTML]FFFFFF${\alpha }_{30}$	\cellcolor[HTML]FFFFFF${\alpha }_{31}$	\cellcolor[HTML]FFFFFF${\alpha }_{32}$	\cellcolor[HTML]FFFFFF${\alpha }_{33}$
NaCl	-88557.123	86351.784	-29910.216	413.046	-7125.857	18640.780	-18244.074
CaCl2	-59216.108	57894.824	-20222.898	217.778	-3770.645	9908.135	-9793.484
\cellcolor[HTML]FFFFFFCon.	\cellcolor[HTML]FFFFFF${\alpha }_{34}$	\cellcolor[HTML]FFFFFF${\beta }_{10}$	\cellcolor[HTML]FFFFFF${\beta }_{11}$	\cellcolor[HTML]FFFFFF${\beta }_{12}$	\cellcolor[HTML]FFFFFF${\beta }_{13}$	\cellcolor[HTML]FFFFFF${\gamma }_1$	\cellcolor[HTML]FFFFFF${\gamma }_2$
NaCl	6335.275	241.57	-980.97	1482.31	-750.98	-0.00134	0.00056
CaCl2	3455.587	307.24	-1259.10	2034.03	-1084.94	-0.00493	0.00231

Equations

$$\mathrm{Ra}$$

$$\mathrm{Sh}$$

$$F_{\mathrm{closed}}\approx F_{\mathrm{open}}{C}_s^2,$$

$$n_{diss.,t}^{\mathrm{C}{\mathrm{O}}_2}={\left(\frac{PV}{zRT}\right)}_i-{\left(\frac{PV}{zRT}\right)}_t,$$

$$\frac{\mathrm{d}{n}_{diss.}^{\mathrm{C}{\mathrm{O}}_2}}{\mathrm{d}t}=A{F}_c,$$

$$F=m^{F.G}{t}_D^{0.5}-b.$$

$$m^{F.G}$$

$$b$$

$$t_{D,\mathrm{NaCl}}^{Q.S}$$

$$t_{D,{\mathrm{NaCl}+\mathrm{CaCl}}_2}^{Q.S}$$

$$F=a^{Q.S}=\left(2.30\times 10^{-8}\right){\mathrm{Ra}}^{1.126}.$$

$$t_{D,\mathrm{NaCl}}^{Sh.}$$

$$t_{D,{\mathrm{NaCl}+\mathrm{CaCl}}_2}^{Sh.}$$

$$\ln x_{{\mathrm{CO}}_2}=\ln \left(y_{{\mathrm{CO}}_2}{\varphi }_{{\mathrm{CO}}_2}P\right)-\frac{{\mu }_{{\mathrm{CO}}_2}^{1\left(0\right)}}{RT}-2{\lambda }_{{\mathrm{CO}}_2-\mathrm{Na}}\left(m_{\mathrm{Na}}+m_{\mathrm{K}}+2m_{\mathrm{Ca}}+2m_{\mathrm{Mg}}\right)-$$

$$-{\zeta }_{{\mathrm{CO}}_2-\mathrm{Na}-\mathrm{Cl}}{m}_{\mathrm{Cl}}\left(m_{\mathrm{Na}}+m_{\mathrm{K}}+m_{\mathrm{Mg}}+m_{\mathrm{Ca}}\right)+0.07m_{SO4},$$

$$\mathrm{Par}\left(T,P\right)=c_1+c_{2}T+\frac{c_3}{T}+c_4{T}^2+\frac{c_5}{\left(630-T\right)}+c_{6}P+c_{7}P\ln T+c_8\frac{P}{T}+c_9\frac{P}{\left(630-T\right)}+c_{10}\frac{P^2}{{\left(630-T\right)}^2}+c_{11}T\ln P$$

$${\varphi }_{{\mathrm{CO}}_2}=c_1+\left[c_2+c_{3}T+\frac{c_4}{T}+\frac{c_5}{T-150}\right]P+\left[c_6+c_{7}T+\frac{c_8}{T}\right]P^2+\frac{\left[c_9+c_{10}T+\frac{c_{11}}{T}\right]}{\ln P}$$

$$\ln \left(\frac{{\mu }_{\mathrm{brine}\mathrm{of}\mathrm{C}{\mathrm{O}}_2}}{{\mu }_{{\mathrm{H}}_2\mathrm{O}}}\right)=Am+B{m}^2+C{m}^3,$$

$$A$$

$$B$$

$$C$$

$$\Delta {\mu }_{mix.}$$

$${\mu }_{mix.}$$

$$\Delta {\mu }_{\mathrm{NaCl}}$$

$$\Delta {\mu }_{\mathrm{CaC}{\mathrm{l}}_2}$$

$${\rho }_b=\frac{1000{\rho }_w+M_{\mathrm{NaCl}}{m}_{\mathrm{NaCl}}{\rho }_w}{1000+A_0{m}_{\mathrm{NaCl}}{\rho }_w+B_0{m}_{\mathrm{NaCl}}^{1.5}{\rho }_w+C_0{m}_{\mathrm{NaCl}}^2{\rho }_w}$$

$${\rho }_b\left(T,P,m\right)$$

$$\ln \frac{P_{\mathrm{ref}}\left(T\right)}{P_c}$$

$${\rho }_{\mathrm{ref}}\left(T,m\right)-{\rho }_0\left(T\right)$$

$$B\left(T,m\right)$$

$$C\left(m\right)$$

$$\frac{{\rho }_0\left(T\right)}{{\rho }_c}$$

$$V_{\mathrm{C}{\mathrm{O}}_2}=37.51-9.585\times 10^{-2}T+8.740\times 10^{-4}{T}^2-5.044\times 10^{-7}{T}^3,$$

Full text

I Introduction

Natural convection emerges when a fluid with higher density overlays a lower density fluid. Convection can significantly enhance the transport of mass, heat, and energy, and is ubiquitous in many natural and industrial areas including weather systems (Doswell III, 1987; Hodges and Thorncroft, 1997), plate tectonics (Elder, 1967; Gurnis, 1988), oceanic currents (Kämpf and Backhaus, 1998), sea water and groundwater aquifers (Van Dam et al., 2009), free air cooling (Fedorov and Viskanta, 1997), electronic devices (Rolin et al., 2014; Hu et al., 2018; Ordonez-Miranda et al., 2016), solar ponds (Suárez et al., 2010), and microfluidic systems (Zehavi et al., 2016). In the context of porous media, natural convection has applications in geothermal energy production (Bejan and Khair, 1985), heat exchangers (Cheng and Minkowycz, 1977), sand separation from oil (Mansour and Ahmed, 2015), and many others. The convection in porous media has recently received renewed attention because of its importance in carbon dioxide (CO2) sequestration in the underground formations as one of the most promising options to stabilize atmospheric CO2 concentrations and hence mitigate the global climate change (IPCC, 2005). Deep saline aquifers have been recognized as the primary target among storage repositories beneath the Earth’s surface, mainly owing to their favorable chemistry, porosity, temperature, pressure, huge capacity, and wide distribution all over the world (Keller et al., 2014). Natural convection controls the dissolution of injected CO2 into underlying brine, which is a key mechanism for the permanent and efficient trapping of CO2. Our goal here is to study the dynamics of convective CO2 dissolution in brine, which not only gives insights into the short- and long-term fate of CO2 injected into the subsurface but may also contribute to the further understanding of convection dynamics in other fields.

During the geological sequestration of CO2 in deep brine-bearing formations (i.e., saline aquifers), buoyant CO2 rises upward until it is confined by an impermeable caprock while spreading laterally beneath. This structurally trapped CO2 is, however, susceptible to leak back to the surface due to the existence of small fractures or faults in the seal. Moreover, in the injection stage the caprock integrity may be compromised by the reservoir overpressurization, which would induce new fractures or cause slip along pre-existing faults in the seal (Szulczewski et al., 2013; Soltanian et al., 2018).

CO2 dissolution in brine (or groundwater) is an important trapping mechanism towards permanent storage of CO2, which would reduce the risk of leakage from imperfect or compromised caprocks (Javaheri et al., 2010; Szulczewski et al., 2013; Martinez and Hesse, 2016; Soltanian et al., 2016a; Newell et al., 2018; Shi et al., 2018; Dai et al., 2018). The diffusion of CO2 molecules into brine initiates the dissolution process (Soltanian et al., 2016b, 2017), which results in a diffusive boundary layer that is more dense than the underlying formation brine and hence prone to density-driven instabilities (Soltanian et al., 2017). Beyond a critical thickness of this boundary layer, fingering instabilities will form. As these instabilities grow, they migrate toward the formation bottom through convection while carrying further dissolved CO2 away from the seal. The underlying lighter fluid at lower CO2 concentrations rises upward at the same time, thereby sharpening the CO2 concentration gradient at the gas-brine interface that accelerates the dissolution rate (Hassanzadeh et al., 2005, 2007; Farajzadeh et al., 2007; Emami-Meybodi et al., 2015). An estimation of the dissolution flux helps to constrain the amount of CO2 that will remain in solution in the subsurface and the amount that is prone to escape.

The dissolution and mixing of CO2 has been well studied in idealized systems of porous convection. The two-phase fluid system that forms following the injection of CO2, where free-phase (gaseous) CO2 overlays the brine-saturated layer with a dynamic and temporally evolving interface in between, is typically simplified into a one-phase system through one of the following approaches. Analog fluid systems are used as an alternative, where the two-phase CO2-brine system is replaced with a single-phase two-layer system composed of water and suitable fluid (e.g., MEG or PPG) that is miscible with water. As such, the analog fluid systems result in non-monotonic density profile for mixture, and cannot represent the partial miscibility (Amooie et al., 2017a, b), density and viscosity profiles, instability strength, and thus the underlying dynamics of actual CO2-brine systems (Backhaus et al., 2011; Hidalgo et al., 2012; Raad and Hassanzadeh, 2015). The two-phase system models can be also simplified by including only the porous layer below the gas-brine interface and representing this interface by (i) a top boundary fixed at CO2 saturation(e.g., Pau et al., 2010; Slim et al., 2013; De Paoli et al., 2017) or (ii) a top boundary of constant flux of CO2 at low rate such that CO2 goes into solution immediately (Soltanian et al., 2016b, 2017; Amooie et al., 2018).

In these systems, dissolution process after the onset of convection can be summarized as follows: First, the fingers grow independently and descend to the bottom of aquifer with small lateral interactions (Slim et al., 2013). However, as time goes on fingers start to merge to their neighboring fingers and create stronger fingers. At this time, the dissolution flux (or mixing rate) increases with time in a *flux-growth * (or dissipation-growth) regime (Slim, 2014; Wen et al., 2018). After the convective regime has well developed, the *quasi-steady state *regime starts (Wen et al., 2018), during which the convective flow brings the fresh brine to the gas-brine interface and the process continues in a *quasi-steady *regime until the brine brought to the gas-brine interface begins to contain dissolved CO2. The latter can happen some time after the descending fingers reach the bottom. The dissolution process then shifts to the *shut-down *regime and continues to reach the maximum dissolution capacity of system at the given pressure, temperature, and salinity (Riaz and Cinar, 2014). For completeness, we gather and summarize the scaling results from previous numerical and experimental studies for these regimes in Table I, which are usually based on the Sherwood number (Sh) and the Rayleigh number (Ra) as:

[TABLE]

where $\Delta\rho$, $k$, $g$, $H$, $\mu$, $D$, $\varphi$, $\Delta c$, and $F$ are respectively density difference, permeability, gravitational acceleration, depth, viscosity, diffusion coefficient, porosity, concentration gradient, and dissolution flux.

What all these idealized systems have in common is that they consider only a single-phase system where free-phase CO2 as well as multiphase processes that could affect the interface dynamics, partial pressure evolution for each phase, CO2 solubility and the associated density increases are absent [13, 15, 16, 23, 24]. In addition, high pressure CO2-brine or CO2-water experiments have been conducted only in blind cells and there is no visual and quantitative data in the literature on CO2-brine convection dynamics in porous media that can represent the geological carbon sequestration in deep, high-pressure saline aquifers. Here, we present, for the first time, the convective dissolution behavior of CO2 in brine via engineering a high-pressure Hele-Shaw cell that allows for both qualitative and quantitative insights into the different flow regimes and the underlying two-phase fluid dynamics. Our visual inspection and scaling analysis of density-driven flow following the dissolution of free-phase CO2 in an underlying porous layer more closely reflect the fate of injected CO2 under a geologic trap, and will advance our fundamental understanding and predictive capabilities of CO2 sequestration.

Saline aquifers can either be open, allowing for the compensation of pressure changes by brine migration, or closed (fault bounded), naturally not allowing for pressure change compensation (Akhbari and Hesse, 2017). In closed aquifers, brine saturation with CO2 (Hewitt et al., 2013; Slim et al., 2013) and likely pressure drop in CO2 free phase (Wen et al., 2018) would limit the CO2 dissolution in brine. Wen et al. (2018) showed numerically that in closed aquifers the dissolution flux in the intermediate, supposedly *quasi-steady *regime is actually not constant, and suggested the following relation between the flux in closed and open aquifers (Wen et al., 2018):

[TABLE]

where $C_{s}$ represents the equilibrium concentration of CO2 at the interface with brine. In this study, confined high-pressure Hele-Shaw cell is designed and implemented in an experimental set-up to investigate the complex behavior of CO2 dissolution in closed aquifers.

Most previous studies, moreover, have considered only pure water (i.e., with no salinity) rather than brine. The latter is typically composed of different salts (Aggelopoulos et al., 2011). Data on brine composition shows the dominance of chloride (Cl-) and sodium (Na+) ions. Among other cations, calcium (Ca2+) is the most frequent one (Wellman et al., 2003; Gaus et al., 2005; Xu et al., 2006; Bacon et al., 2009; Azin et al., 2013; Mohamed et al., 2013; Wang et al., 2016; Vu et al., 2017; Shi et al., 2018). Transition time between the dissolution regimes will control the total amount of dissolved CO2 with time. One of our goal here is to provide scaling relations for transition times between dissolution regimes and for dissolution fluxes in CO2 brine systems that contain NaCl and NaCl+CaCl2 mixtures. These relations are helpful for examining system behavior in more realistic conditions where a mixture of salts constitutes the formation brine.

The present study is organized as follows. Section II presents the methodology and the description of experimental set-up and details of the implemented tests. Results are provided and analyzed in section III. In the first part of this section, the dynamic of dissolution process is examined through the visual data. In the second part, quantitative data (pressure data as well as dissolution flux) is studied, and scaling relations are obtained for transition times between dissolution regimes and for dissolution flux. Finally, we draw the conclusions in section IV, followed by Appendix A–E that provide the details for thermodynamic and physical properties of our fluid system.

II Methodology

We conducted series of eight experiments in high-pressure conditions at different medium permeability and salinity values. The details of the experiments conducted are given in Table II. To understand the effect of brine composition we reduced the NaCl amount to 80$\%$ of its initial value, and for the remaining 20$\%$ we added CaCl2. For example, in the corresponding solution of 1 Molal (M hereafter) NaCl, we used 0.8 M NaCl and 0.1 M CaCl2. This way there is similar concentration of Cl- in both solutions, and for each divalent cation of Ca2+ there are two monovalent cations of Na+ in the corresponding solution.

Experiments have been performed in constant volume of a Hele-Shaw cell that we have designed. The cell is of the internal dimensions of 36$\times$30$\times$2.5 cm (internal volume of 2700 cm3) with aluminum frame covered with a 5 cm thickness plexiglass enabling the visual examination (the set-up shown in Figure1). To prevent plexiglass expansion under high-pressure conditions, we attached a steel frame to the main frame. We packed $\sim$ 23 cm from of total height of 30 cm with glass beads. Permeability range is 400–550 Darcy (D), and porosity is $\sim$ 0.36 for all cases. We used a wider opening (2.5 cm) in our cell, which allows for adding such porous structure into the cell as an improvement over previous studies. For instance, the onset of convection in bulk fluid experiments happens so instantaneously that usually there is not sufficient data in the early, diffusion-dominated regime. In the presence of porous structure (system with macroscale tortuosity and lower permeability and porosity than bulk fluid), conversely, we can obtain pressure data from diffusion-dominated regime and hence calculate diffusion coefficient from experimental data. Further, our new device has negligible boundary effects due to its size, as compared to Hele-Shaw cells used previously.

In our experiments, the porous part of cell is saturated with desired brine solution; we introduce brine into the cell by continuous mixing and gradually adding NaCl or mixture of NaCl and CaCl2 into water. Bromocresol green is used as a pH indicator with 0.02% weight fraction in brine. CO2 is injected to the remaining 7 cm in the empty upper section—with no porous structure—to reach the desired pressure. This packing configuration helps to eliminate the capillarity effects as there is no longer a capillary transition zone between the initial brine- and CO2-saturated layers. In practice, The cell is placed in an oven and vacuum pump is used to saturate it with brine and the pH indicator solution. Then, the transfer vessel cylinder is loaded with CO2 gas phase. The CO2-containing cell is pressurized, and all equipment in the oven is kept under the 50∘C for at least 12 hours. Before inserting gas into the cell, a pump vacuumizes the cell. A valve between gas cylinder and the cell is opened gradually to reduce possible disturbances. Once we obtain a desired pressure, we close the cylinder valve and open the pressure gauge valve. Initial pressures for the tests are in the range of 502.6–535.3 psi. Because our set-up has constant volume in the isothermal condition ($\sim$ 50 C), the CO2 dissolution into the brine results in some pressure reduction, the continuous monitoring of which can provide a quantitative measure of dissolution. Pressure changes are therefore measured with a high-accuracy pressure gauge ($\pm 0.1$ psi resolution). The pressure monitoring and the continuous capturing of visual snapshots with a digital camera are performed until pressure reaches to an almost constant value.

We use non-dimensional time as $t_{D}={tD}/{H^{2}}$ throughout this study, where $t$, $D$, and $H$ are time, diffusion coefficient and the height of system, respectively. Based on this non-dimensional time, we present an example picture of the dynamic evolution of pressure resulting from our experiments in **Figure 2 **.

As mentioned above, we calculate the diffusion coefficient values from the pressure data in the diffusion-dominated regime in our experiments, and use these values in this study. By matching the experimental data with simple simulation models following our previous work (Mahmoodpour and Rostami, 2017; Mahmoodpour et al., 2018), we were able to obtain diffusion coefficients. Unlike in previous experiments with blind cells where the gas-brine interface is unknown (e.g., Mojtaba et al., 2014), the visual data from our experiments allow us to obtain the volume of CO2 gas phase ($V$). We use the equation below to calculate the amount of dissolved CO2 (Duan and Sun, 2003):

[TABLE]

where $n_{\mathrm{diss.},t}^{\mathrm{CO}_{2}}$ shows the amount of dissolved CO2 and i stands for the initial time condition. The mole fraction of brine in the gas phase is assumed negligible (see (Duan and Sun, 2003)).

The accurate detection of different dissolution regimes is not feasible by pressure data only. Therefore, the dissolution flux $F_{c}$ is used instead, which is calculated using the values of dissolved CO2 in brine as follows:

[TABLE]

where ${A}$ is the surface area of the gas-brine interface (here 0.025$\times$0.36 m2).

To detect the *quasi-steady-similar * regime in our closed system at the constant volume condition, we follow the analytical approach by Wen et al. (2018), with a difference that we use the $C_{s,0}^{2}$ (initial concentration of CO2 at interface) in our conversion to make it non-dimensional as $\frac{F_{c}}{C_{s}^{2}}\times C_{s,0}^{2}\approx F$ that can be used for obtaining scaling relations. We will discuss this modification later with details. Since the temperature and salinity are constant, the available correlations for the solubility of CO2 in brine can be used to calculate the Cs values at any given pressure. Here, we use the solubility models by Duan and Sun (2003) and Duan et al. (2006) to calculate the CO2 equilibrium concentration at the gas-brine interface. Importantly, these models show that an increase in salinity will result in a reduction of solubility (i.e., equilibrium concentration).

Figure 3 shows the resulting dissolution flux for the presented pressure curve in Figure 2. It should be noted that, the fluctuating data of dissolution flux are converted to produce a smooth curve via moving average method. Based on the dissolution flux dynamics, the transition between dissolution regimes is detectable. The lines illustrated in the figures show the onsets of convection, quasi-steady (-similar), and shut-down regimes, respectively, from left to right.

While a pressure gauge with a high accuracy of $\sim$ 0.1 psi was used in this study, pressure changes due to CO2 dissolution over a small time period on the order of seconds are not detectable. We analyzed the pressure data based on eight-minute intervals.

III Results and Discussion

Below we present our results for the density-driven flow of CO2 in brine solutions with dissolved NaCl and NaCl+CaCl2. Specifically, we analyze the dissolution flux and the critical times for different flow regimes. Scaling relations are also presented. To verify the presented scaling relations against experimental observations, the individual transition times and dissolution rates are further compared with the predicted values.

III.1 Dynamics of Dissolution

At early time, the dissolution of CO2 into brine is a diffusion-dominated process. The mechanism behind the diffusive regime is well-documented in the literature; the dissolution flux exhibits the classical Fickian scaling of $F\propto t^{-0.5}$ (Szulczewski et al., 2013; Soltanian et al., 2016b, 2017). We present in Figure 4 the visual results for our Case 4 (introduced in Table II), where the brine has 1 *M *of NaCl and the permeability is 400 D. Figure 4a shows the first diffusive regime, where the diffusive layer can be seen as an almost piston-like descending layer near the gas-brine interface.

The more dense CO2-rich brine overlying the lighter fresh brine leads to a gravitationally unstable density stratification. Diffusion has damping effect on these instabilities. As per Szulczewski et al. (2013), in conditions with Ra greater than some critical values (e.g., $Ra>55$) instabilities grow with time and result in the development of convective fingers (see Figure 4b). At this time, although the convective fingers are clear, their strength is not sufficient to make substantial changes to the pressure. In early times, convective fingers grow independently. As time passes and convective fingers grow enough, there will be strong interactions between fingers.

The mechanisms that control the interactions between gravitational fingers following the dissolution of CO2 in our brine-saturated porous medium are observed as follow. 1) Side merging depicted by a square box in Figure 4c: the downward motion of a dominant finger as well as the upward migration of fresh brine from bottom create strong circulating velocity field toward the center of finger, which in turn drives smaller fingers toward the dominant finger. Eventually, smaller fingers merge into the dominant finger, thus increasing its strength. 2) Root zipping depicted by (solid) rectangular box in Figure 4e and Figure 4f: two fingers merge from their roots. 3) Tip splitting depicted by an oval in Figure 4f: the tip of a convective finger is flattened by the upward motion of brine and is likely to split into different portions. 4) *Necking *depicted by a circle in Figure 4j: in some part of a finger the width decreases and its feeding from the upper part decreases. 5) Trailing lobe detachment depicted by a dashed oval in Figure 4k: a portion of the convective finger is separated from the root. 6) *Reinitiation *depicted by dashed circles in Figure 4j and Figure 4k: new small-scale fingers, known as protoplumes, are initiated between the dominant descending fingers in intermediate times during dissolution with the generation of concentration gradients below the interface due to upwelling flow of fresh water. 7) *Fading fingers which partly is clear in *the dashed rectangle in Figure 4d: these fingers are created in early times and grow in a similar way as the surviving fingers. At some point, the growth of fading fingers cease, either when they intersect dominant plumes or through diffusive smearing, and they will act as sourcing pool for feeding the dominant fingers. Remarkable similar interactions have been reported in prior modeling work on viscous fingering and diffusive-convective dissolution (Tan and Homsy, 1988; Zimmerman and Homsy, 1991, 1992; Ghesmat and Azaiez, 2008; Ranganathan et al., 2012; Amooie et al., 2017b, a; Sabet et al., 2018; Amooie et al., 2018).

Under the above mechanisms the dominant convective fingers grow larger and stronger to reach a point at which the pressure decay data show a rapid decline (in comparison to the diffusion-dominated regime). This point is selected as the onset of convection, beyond which the pressure decreases at a greater rate in comparison to the diffusive regime. This time is an important factor from the operational stand-point since there is pressure decline below the caprock. which implies a smaller magnitude of stress and thus a smaller risk of induced leakage.

After the onset of convection there are two main mechanisms that enhance the dissolution flux. First, as the fingers of CO2-rich brine sink down, fresh brine occupies their place due to mass conservation. Therefore, the concentration gradient at gas-brine interface becomes steeper and mass transfer increases. Second, when fingers sink down, the surface area between sinking fingers and rising fresh brine increases while plume stretching simultaneously steepens the concentration gradients in the direction perpendicular to the finger, the combined effect of which leads to mass transfer enhancement and dissolution across the formation (Borgne et al., 2014). This forms the flux-growth regime that follows the onset of convection, where the dominant fingers (or megaplumes) start to appear and will further develop, thereby enhancing the dissolution flux (see Figures 4f, 4g, 4h, and 4i).

To better elucidate the behavior of the convective fingers resulting in the flux enhancement, a narrow slice is selected at the depth of 5 cm below the gas-brine interface and the formation of dissolved CO${2}$ patterns over time is captured by a representative spatiotemporal diagrams of the dense fingers (Figure 5). Except for Case 6, all cases follow a similar trend: the number of convective fingers and their speed increase with Ra. These characteristics (fingering number and their speed) are stronger for NaCl solutions in comparison to their corresponding NaCl+CaCl2 mixtures. This is fundamentally because of the smaller diffusion coefficient for the former (NaCl) than the latter, which would result in larger Ra and thus more pronounced fingering (note that density difference and viscosity profiles for the corresponding brines of similar salinity are indeed similar to each other). The characteristic features of coarsening dynamics that have been identified in previous numerical studies and experiments on analog fluids are interestingly captured in our CO2-brine experiments at high pressure as well. As an example, the Case 3 of Figure 5 shows that when convective fingers sink down, they create low-resistance pathways where new high-concentration fluid preferentially follows the long-lived descending plumes towards the bottom of system (signified by an oval shape symbol in the figure). As indicated by circle symbols in the Case 4 of Figure 5, the merging between fingers is more predominant early on when small-scale instabilities emerge. Eventually, the continuous merging between neighboring plumes leads to the development of few larger-scale coherent structures that serve as the conduits for traveling plumes (white color circle symbol). The rectangular box in the Case 4 and 7 represent those weak convective fingers that fade while separating from a dominant one and diffusing into the main stream.

Traditionally, when there is no gas cap and pressure change subject to dissolution, the flux grows to a maximum (owing to the mechanisms explained above), beyond which merging and shielding between adjacent elongated fingers combined with diffusive spreading weaken the concentration gradients in the boundary layer. Hence the dissolution flux stops to grow, and in fact it will decrease transiently until a quasi-steady constant-flux regime develops as a result of a balance between the reinitiation of protoplumes and the coarsening of the existing plumes followed by their subsumption (Amooie et al., 2018). Our experiments, however, reveal that such quasi-steady constant-flux regime actually never happens in a two-phase closed system, as suggested by Figure 3. The pattern coarsening, as the mechanism for flux reduction, dominates the protoplume reinitiation that is the mechanism for flux enhancement (denoted by dashed circles in Figures 4j and 4k). This is because of the predominant decrease of pressure following the onset of convection, which in turn decreases the solubility of CO2 in brine at the interface and accordingly the dissolution-induced density change, i.e., the driving force for buoyancy-driven convection. The final outcome of this nonlinear dynamics is the decreasing trend, rather than a constant average value, of flux in the intermediate regime before the onset of final shut-down, as predicted theoretically by Wen et al. (2018).

Despite the decreasing trend of flux in this regime, a quasi-steady regime, however, is attainable in closed systems with the introduced compensation of dissolution flux as $\frac{F_{c}}{C_{s}^{2}}\times C_{s,0}^{2}\approx F$—shown in Figure 3. This is a modification to the scaling suggested by Wen et al. (2018) as ${F_{c}}/{C_{s}^{2}}$. The results of experiments are analyzed based on the new correction to obtain transition times between dissolution regimes and the associated scaling relations.

Figure 6 shows the fingering patterns at the time of transition to the quasi-steady-similar regime for all cases. This regime persists until the upwelling flow starts to carry CO2-rich fluid towards the gas-brine interface. Next, the concentration gradient at the interface and consequently the dissolution flux decays rapidly with time, resulting in the shut-down of convective dissolution process. The transition to the convection shut-down regime is presented in Figure 4n, after which the position of megaplumes is almost constant and the region between adjacent megaplumes is slowly saturated mainly by diffusion mechanism (Figure 4o).

Our findings suggest that the onset of shut-down regime occurs after twice the time required for the first fingers to reach the bottom. The results show that the first contact, e.g., for Case 4 occurs at almost 480 min after the start of experiment (Figure 4h), but the shut-down regime starts at approximately $1160>2\times 480$ min. Further analysis of the visual data explains this behavior. As displayed in Figures 4j and 4l with the dashed rectangular boxes, when the convective fingers reach the impermeable bottom of system they initially propagate laterally therein before rising upward with the upwelling flow. The traveling speed of rising wavefront following the impact with bottom boundary is slower than the sinking plumes, partly due to a greater lateral extent (surface area) for the rising fluid (i.e., a mixture of brine and dissolved CO2) as compared to each sinking finger. Therefore, the onset time of shut-down regime is larger than the twice the time required for the first impact with bottom boundary. Note that blind-cell experiments could not observe this. It should be also noted that the transition times between dissolution regimes are obtained here through quantitative data and will be discussed below with details.

III.2 Scaling Analyses

Pressure data obtained from our tests are presented in Figure 7. During the diffusion-dominated regime pressure (besides temperature) is almost equal for all cases. Therefore, salinity controls the diffusion coefficient and consequently the dissolution rate during the first diffusion-dominated regime (see Table II). NaCl+CaCl2-containing brines show a slightly higher diffusion coefficient than the corresponding NaCl brines, most likely due to the lower total molality for NaCl and CaCl2 mixture as well as the different molecular interactions in different salt solutions (see Table II). In light of this observation, the dissolution rate for NaCl+CaCl2 brine is higher in the first diffusion-dominated regime. The increase in salinity results in the reduction of diffusion coefficient, which leads to the lower dissolution rates in the diffusion-dominated regime, and thus lower pressure decline, for cases with higher salinity (Figure 7).

To obtain scaling relations for dissolution flux and the transition times between dissolution regimes, we converted the pressure data to dissolution flux with the help of visual data for the dynamic volume of gas phase. We calculated the modified dissolution flux based on Eq. 3 (Wen et al., 2018) for the quasi-steady-similar regime to investigate the possible scaling relations with this character. Note that such modified dissolution flux in the quasi-steady-similar regime is independent of time. Results are reported in Table III.

We find that ${F_{c}}/{C_{s}^{2}}$ surprisingly does not show any specific trend with Ra, unlike the predictions of Wen et al. (2018). To better understand the reason behind this, let us consider two situations: i) Ra is increased through an increase in permeability (e.g., Case 1 versus Case 2). In this situation the dissolution flux (${F_{c}}$) increases with Ra at a given time while there are no considerable changes in the CO2 solubility values ($C_{s}$); therefore, the overall value of ${F_{c}}/{C_{s}^{2}}$ is found to increase with Ra in this situation (e.g., from 5 in Case 2 to 7 in Case 1). ii) Ra is increased through a reduction in salinity (e.g., Case 1 versus Case 3). Whereas the ${F_{c}}$ increases in this situation, the $C_{s}$ is increased as well (see Section II), and the overall value of ${F_{c}}/{C_{s}^{2}}$ actually decreases (from 7 in Case 1 to 4.1 in Case 3). Based on this explanation and opposing situations, we find that how Ra changes specifically (via salinity or permeability) is important for the dynamical behavior of ${F_{c}}/{C_{s}^{2}}$, and there is not a general relation between the latter and Ra.

Given the definition of Sherwood number as $Sh\!=\!\frac{F_{c}H}{D\varphi\Delta C/H}\!=\!\frac{F_{c}H}{D\varphi c_{b,\mathrm{max}}^{s}x_{\mathrm{CO}_{2},\mathrm{max}}^{s}/H}$ (Amooie et al., 2018), Sh values show decreasing trend during the quasi-steady-similar regime in a closed system. This is indeed because of the decreasing dynamics of ${F_{c}}$ in this regime as shown previously and the other parameters assumed as constant (see Figure 3). Further, $c_{b,\mathrm{max}}^{s}$ (constant maximum density of brine) and $x_{\mathrm{CO}_{2},\mathrm{max}}^{s}$ (constant maximum solubility fraction) in the formulation of Sh do change with salinity, and thus the above-mentioned problem with ${F_{c}}/{C_{s}^{2}}$ scaling with Ra applies to Sh scaling as well.

To overcome these physical challenges in obtaining rigorous scaling relations, we propose a further adjustment as $\frac{F_{c}}{C_{s}^{2}}\times C_{s,0}^{2}\approx F$, where $C_{s,0}^{2}$ is the initial concentration of CO2 at interface and itself sensitive to salinity, such that the dynamic changes in $C_{s}^{2}$ will be compensated with $C_{s,0}^{2}$. The results are presented in Figure 8 for all Cases.

The analysis of dissolution flux dynamics shows that for both the solutions of NaCl and NaCl+CaCl2, the dissolution flux increases with time as $t_{D}^{0.5}$ in the flux-growth regime. This is consistent with the findings of Newell et al. (2018). Higher Rayleigh numbers result in a greater number of convective fingers with more strength. Therefore, we hypothesize that the dissolution flux at a given equal time in the flux-growth regime should be greater than that for a lower Ra case, i.e., the pre-factor for the $t_{D}^{0.5}$ scaling should be positively correlated with Ra. Furthermore, as shown in Figure 5 and in supplementary videos (See Supplemental Material at [URL will be inserted by publisher] for the videos of laboratory experiments on convective dissolution., ), the number of convective fingers and their strength (quantified with sinking speed) are higher for brine with NaCl than for brine with NaCl+CaCl2 mixtures. Dissolution flux data can quantify these observations. The relation with the following formula are used to fit the dissolution flux data (with continuous curve):

[TABLE]

Values of the $m^{F.G}$ and $b$ are reported for each experiment denoted by $m_{fit.}^{F.G}$ and $b_{fit.}^{F.G}$ in Table III. Detailed analysis of these factors prove that they interestingly scale linearly with Ra (Eqs. 7 and 8). The estimated values based on these equations are further reported in Table III as $m_{est.}^{F.G}$ and $b_{est.}^{F.G}$ .

[TABLE]

The modified dissolution flux continues to increase until it plateaus in the quasi-steady period. This time is detected as the onset of quasi-steady regime ($t_{D,fit.}^{Q.S}$), and is reported in Table III. Results reveal that this onset occurs earlier for cases with higher Ra in both mixture types. This is because the interactions between convective fingers are more intense at higher Ra, and therefore the megaplumes form earlier during the experiments. Based on the same reasoning, NaCl solutions show earlier onset for the quasi-steady regime in comparison with their corresponding NaCl+CaCl2 solutions. Detailed analysis of the onset of quasi-steady regime show that they correlate through the following equations with Ra in two mixtures:

[TABLE]

Estimated times for onset of quasi-steady regime ($t_{D,est.}^{Q.S}$) with the proposed scaling relations are presented in Table III, which show good agreements with the fitted values. As discussed in Figure 4, protoplume reinitiation helps to maintain the dissolution flux during quasi-steady regime. Cases with higher Ra show higher (modified) dissolution flux during this regime. These constant dissolution rates (which are presented with $a_{fit.}^{Q.S}$ ) are fitted with the Eq. 11 for all cases and the resulting estimations are reported with the symbol of $a_{est.}^{Q.S}$ .

[TABLE]

The almost linear scaling of compensated dissolution flux with Ra obtained from our novel experiments is somewhat striking given the debate that exists for the flux scaling (linear versus sublinear) in the quasi-steady regime (see Table I).

Dissolution process continues in the quasi-steady regime until the wavefront of CO2-rich brine from the bottom reaches the interface. It is expected to find an earlier onset of shut-down regime for higher Ra cases just because of its control on the traveling speed of convective fingers. We observe that our experimental results meet this expectation. Another key feature we find is that NaCl+CaCl2 solutions show longer duration of quasi-steady regime as compared to the corresponding NaCl solutions (Figure 8). This can be justified by our previous explanation as to why the onset of shut-down takes more than two times the time of the first impact with bottom boundary; lower Ra for the NaCl+CaCl2 solutions compared to their corresponding NaCl brine would result in even longer travel times for the rising wavefront. Hence, the onset of shut-down occurs much later for the NaCl+CaCl2 than for an equivalent NaCl solution, accordingly with longer duration of the quasi-steady regime. The analysis of the onset of shut-down times ($t_{D,fit.}^{Sh.}$) shows that following relations describe well the transition times. The corresponding estimated values ($t_{D,est.}^{Sh.}$) with these relations are reported in Table III and show close agreement with the fitted values.

[TABLE]

During the shut-down regime, the CO2 concentration gradient at gas-brine interface constantly decreases while reducing the dissolution flux. Whereas the previous theoretical and numerical studies reported a $F\propto t^{-2}$ characteristic behavior during this regime (Slim, 2014; De Paoli et al., 2017; Wen et al., 2018), recent experiments of Newell et al. (2018) suggested a $F\propto t^{-1.75}$ scaling during the shut-down regime. Our experiments remarkably confirmed the temporal scaling behavior of the theoretical and numerical studies, where $F=a^{S.D}t_{D}^{-2}$ describes our results (see the continuous line during this regime in Figure 8). Pre-factors of $a^{S.D}$ are also reported in Table III. These values, however, do not show a distinct trend with Ra, and it seems they are a function of diffusion coefficients and depend on the history of prior regimes. Cases with higher salinity show larger values for pre-factor, because of their lower concentration of dissolved CO2 as a result of prior dissolution events. The NaCl+CaCl2 mixtures show higher pre-factors in comparison to their corresponding NaCl solutions because of larger diffusion coefficients.

IV Discussion and Conclusions

In this study, high-pressure laboratory Hele-Shaw experiments are performed to examine the convective dissolution of actual CO2 in a closed porous media system saturated with brines that contain NaCl and NaCl+CaCl2 mixtures at different molality. The common concerns of previous studies that use analog fluid systems to emulate the CO2-brine behavior, such as non-monotonic density profile and miscibility of system, are resolved in our work. Further, in comparison to the experiments performed in blind cells, visual data are available in this study enabling us to examine the dynamics of convective fingers during dissolution. Performing the experiments at high pressure not only more closely reflects the subsurface conditions but also makes it possible to obtain diffusion coefficients and dissolution fluxes through the measurable pressure changes during the dissolution process. Our new engineering provides an unprecedented opportunity to bridge the gaps among the visual low-pressure and high-pressure blind cell bulk-fluid experiments and the analog fluid experiments in porous Hele-Shaw cells at atmospheric conditions. Juxtaposing the qualitative and quantitative data in this study further decreases the possible errors during data analysis.

Qualitative data (image) analyses reveal that in the early times, the dissolution process is diffusion-dominated. The dissolution of CO2 into underlying brine creates a flat diffusive boundary layer under the gas-brine interface. The mixture of CO2 and brine overlying the lower density fresh brine leads to a gravitationally unstable stratification. When the diffusive boundary layer is thickened enough, instabilities initiate and grow as convective fingers. Following their appearance, convective fingers grow independently, early on, but as time goes on they interact more with each other and create stronger fingers. Seven mechanisms for fingering interactions are identified in this study as follows: side merging; root zipping; tip splitting; necking; trailing lobe detachment; protoplume reinitiation; fading fingers. We found that the number of convective fingers, their speed and interactions are higher in NaCl solutions as compared to the corresponding NaCl+CaCl2 solutions. In each type of solutions, the speed of fingers and their interactions are positively correlated with Ra. Increasing salinity results in a decrease of convective flow strength. After some period of time, convective fingers well develop into larger-scale coherent structures and quasi-steady-similar regime starts where new smaller scale fingers emerge at the gas-brine interface helping the maintenance of dissolution flux. However, the considerable reduction in gas pressure following CO2 dissolution, especially after the onset of convection, results in a reduction in solubility and concentration at the gas-brine interface. The flux-enhancing protoplume reinitiation is thus reduced, such that it is no longer sufficient to balance the flux-reducing coarsening mechanisms. Final outcome is a progressively decreasing trend in flux within the intermediate regime before the onset of shut-down regime. The latter begins when the brine brought to the interface starts to contain dissolved CO2 in late times, following which the position of megaplumes is almost constant and the region in between is slowly saturated mainly through diffusion. The onset of shut-down regime occurs earlier for NaCl solutions.

Our experiments are conducted in constant volume condition. Consequently, the dissolution of CO2 into brine is accompanied by a reduction of the pressure of gas phase. The analysis of pressure data provides a unique quantitative measure of the dynamic process. To detect the transition between dissolution regimes and reveal the scaling behavior of dissolution, the dissolution fluxes are calculated from the pressure curves with the help of visual data for gas volume. Further, we introduced a modification to dissolution flux in order to yield insights into the scaling of quasi-steady regime (similar to open systems) as a function of Ra. With the conversion of $\frac{F_{c}}{C_{s}^{2}}\times C_{s,0}^{2}\approx F$, dissolution fluxes in the closed system appear to become constant during the quasi-steady regime. Using this definition of flux, four dissolution regimes are identified as follows: diffusion-dominated; flux-growth; quasi-steady and shut-down regimes. Scaling relations based on the dimensionless numbers are introduced for the transition times between these regimes. All transition times for the quasi-steady and shut-down regimes for different salt types exhibit an inverse proportionality with Ra. However, the comparison between results reveals that NaCl solutions show earlier onsets of quasi-steady and shut-down regimes as compared with those from the corresponding NaCl+CaCl2 solutions. The difference between the onset times of shut-down in the two types of brine is larger than those of the quasi-steady regime. Therefore, NaCl+CaCl2 systems experience a longer period in the quasi-steady regime. Dissolution flux in all cases follows the $F\propto t_{D}^{0.5}$ and $F\propto t_{D}^{-2}$ scaling behavior for the flux-growth and shut-down regimes, respectively. During the flux-growth regime, the pre-factors for $F\propto t_{D}^{0.5}$ for NaCl solutions are higher than those for NaCl+CaCl2 solutions, and both linearly scale with Ra. Further, NaCl solutions show higher dissolution fluxes than NaCl+CaCl2 mixtures during the quasi-steady regime, with the compensated dissolution fluxes almost linearly scaling with Ra irrespective of salt types. For the shut-down regime, the dissolution rates for NaCl+CaCl2 mixtures are also higher.

Findings of this study advance our knowledge of CO2 dissolution into formation brine following geological sequestration in deep saline aquifers. Scaling relations are introduced based on the dimensionless numbers, which enable us to use them in similar studies and possibly field experiments. Since CaCl2 is the most common salt after NaCl in the majority of aquifers, this work helps to obtain a more realistic estimation for the short- and long-term fate and transport of CO2—particularly in terms of the dissolution flux and transition times between critical regimes—following its storage in the subsurface.

Appendix A Solubility of CO2 in brine

Solubility of CO2 in brine can be obtained from Duan and Sun (2003); Duan et al. (2006):

[TABLE]

where T in K and pressure in bar. $\mu_{{\mathrm{CO}_{2}}}^{1\left(0\right)}$ , $\lambda_{{\mathrm{CO}_{2}}-\mathrm{Na}}$ and $\zeta_{{\mathrm{CO}_{2}}-\mathrm{Na}-\mathrm{Cl}}m_{\mathrm{Cl}}$ are standard chemical potential, interaction parameter between CO2 and Na+ and interaction parameter between CO2, Na+ and Cl- respectively. Following equation is used to calculate these parameters:

[TABLE]

The fugacity coefficient ($\phi_{{\mathrm{CO}_{2}}}$) is calculated from the following non-iterative equation:

[TABLE]

Constants of c1- c15 are presented in Table IV.

Appendix B Viscosity of NaCl solution (Mao and Duan, 2009)

[TABLE]

where A, B and C are functions of temperature (in K):

[TABLE]

c1- c8 are obtained from the experimental data and available in Table V.

Appendix C Viscosity of CaCl2+ NaCl solution

To obtain the viscosity of the NaCl+ CaCl2 solution, the simple mixing rule is used (Zhang et al., 1997):

[TABLE]

where $c_{1}$–$c_{10}$ are calculated from the experimental data on the specific concentration of the salts. In the required concentrations, parameters are interpolated. For 1 M. concentration of brine, these values are presented in Table V.

Appendix D Brine density

NaCl solution density (Potter and Brown, 1977):

[TABLE]

To obtain the NaCl+ CaCl2 solution density, the relative of the NaCl+ CaCl2 mixture to the NaCl solution is calculated by Al Ghafri et al. model (Al Ghafri et al., 2012):

[TABLE]

where Tc and Pc are critical values and $\varphi=1-{T}/{T_{c}}$, and

[TABLE]

where ${c_{1}}=1.992741,{c_{2}}=1.099653,{c_{3}}=-0.510839,{c_{4}}=-1.754935,{c_{5}}=-45.517035,{c_{6}}=674694.45,{\beta_{00}}=-1622.40,{\beta_{01}}=9383.80,{\beta_{02}}=-14893.80,{\beta_{03}}=7309.10,{\gamma_{0}}=0.11725,{\sigma_{1}}=-7.859518,{\sigma_{2}}=1.844083,{\sigma_{3}}=-11.786650,{\sigma_{4}}=22.680741,{\sigma_{5}}=-15.961872,$ and $\sigma_{6}=1.801225$. Required parameters for NaCl and CaCl2 are presented in Table VI.

Appendix E Partial molar volume of CO2 in water (Garcia, 2001)

[TABLE]

where ${T}$ is in degrees Celsius and $V_{\mathrm{CO_{2}}}$ is in cm3/mole.