Quantifying the dynamic density driven convection in high permeability packed beds

Quantifying the dynamic density driven convection in high permeability packed beds

Magnetic Resonance Imaging 39 (2017) 168–174 Contents lists available at ScienceDirect Magnetic Resonance Imaging journal homepage: www.mrijournal.c...

2MB Sizes 0 Downloads 10 Views

Magnetic Resonance Imaging 39 (2017) 168–174

Contents lists available at ScienceDirect

Magnetic Resonance Imaging journal homepage: www.mrijournal.com

Original contribution

Quantifying the dynamic density driven convection in high permeability packed beds Ying Teng a, Lanlan Jiang a,b,⁎, Yingting Fan a, Yu Liu a, Dayong Wang a, Abuliti Abudula c, Yongchen Song a a b c

Key Laboratory of Ocean Energy Utilization and Energy Conservation of Ministry of Education, Dalian University of Technology, 116024 Dalian, China Research Institute of Innovative Technology for the Earth, Kizugawa City, Kyoto 619-0292, Japan North Japan New Energy Research Center, Hirosaki University, Aomori City 0300831, Japan

a r t i c l e

i n f o

Article history: Received 15 February 2017 Received in revised form 9 March 2017 Accepted 11 March 2017 Available online xxxx Keywords: Density driven convection Finger growth regimes Onset time Mixing time Rayleigh number

a b s t r a c t The density driven convection phenomenon is expected to have a significant and positive role in CO2 geological storage capacity and safety. The onset and development of density-driven convective on the core scale is critical to understand the mass transfer mechanism. In this paper, laboratory experiments were conducted to investigate the density-driven convective in a vertical tube. The deuterium oxide (D2O)/manganese chloride (MnCl2) water solution in water or brine were as an analog for CO2-rich brine in original brine. Experiments are repeated with variations in permeability to vary the characteristic Rayleigh number. Based on the MRI technology, the intensity images showed the interface clearly, reflecting the transition from diffusion to convective. With the echo-multislice pulse sequence method, the intensity images can be obtained as 2 min 8 s. For the denser fluid pairs, fingers appeared, propagated, coalesced and multi-fingers formed. The finger growth rate of the convective was visualized as three distinct periods: rising, stable and declining. Detailed information regarding the wave number, wave length, onset time and mixing time as functions of Rayleigh number are developed. © 2017 Elsevier Inc. All rights reserved.

1. Introduction Carbon capture and sequestration (CCS) is expected to play a major role in meeting CO2 emission reduction targets [1]. The injected CO2 formed a thin layer of free phase below the caprock; and then the above CO2 dissolve in the brine by molecular diffusion creates a dense boundary layer of CO2-rich brine in comparison to CO2-poor brine densities [2]. The density of the boundary layer at the interface is 0.1–1% denser than the original brine, depending the on pressure, temperature and salinity [3]. Due to the gravitational instability, a pattern of the downward flow of CO2-rich brine that looks like fingers will appear, referring to as convective fingers [4]. In contrast to molecular diffusive dissolution, the so-called density-driven natural convection significantly increases the rate of CO2 dissolution, reduces time for dissolution trapping and enhances the safety of the storage in the saline aquifer [5]. Density-driven natural convection has received increasing attention both from theoretical [6,7] and experimental [8,9] point of view. The earliest report [10] investigated density-driven natural convection effect on CO2 dissolution, and then the convection stability criteria for ⁎ Corresponding author at: Key Laboratory of Ocean Energy Utilization and Energy Conservation of Ministry of Education, Dalian University of Technology, 116024 Dalian, China. E-mail address: [email protected] (L. Jiang).

http://dx.doi.org/10.1016/j.mri.2017.03.004 0730-725X/© 2017 Elsevier Inc. All rights reserved.

vertical flow that subjected to both thermal and density gradients [11] and for onset time in the anisotropic porous media [12]. Besides of the numerical study, several laboratory scale experiments were conducted to analyze convective dissolution by using Hele-Shaw cell [8]; [13–15], and the convective mixing process with optical techniques, such as a Schlieren technique [16] and or holographic interferometry [17]. These technology of Hele-Shaw cell provides visualization for convective fingers propagate, however, only mathematically analogous to flow in porous media [18] and limited at ambient conditions [19]. Then the amount of CO2 dissolved into solution, bulk dissolution rate and the onset of the density-driven convection were also investigate at high pressure by using mercury-free DBR pressure-volume temperature (PVT) [20], stainless-steel cylindrical blind cell [2,21] by monitoring pressure decay. The results showed the mass transfer rate was usually larger than expected from a diffusion process alone. Those experiments cannot provide visual observation of occurrence of convective and the instabilities because of using blind cell. The dimensionless parameter Rayleigh number (Ra), nominally expresses the ratio of convective and diffusive transport in the natural convection process [22]. The density-driven natural convection is expected for Rayleigh number above 4π2 in porous media [23]. While the mass transfer or convective flux [19], the onset time of convective currents [24], the critical wave length and the wave numbers of the finger [25] are all as functions of the Rayleigh number in the range of

Y. Teng et al. / Magnetic Resonance Imaging 39 (2017) 168–174

169

Fig. 1. Schematic diagram of the experimental setup. Left shows the MRI experimental machine with the connector tube. Right shows the convection for one fluid pair in BZ-06. Denser fluid was injected from the bottom and indicted with blue color. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

5000–30,000 [26], with similar conditions in real saline aquifers, due to the different properties of porous media [27]. The visual observation for investigating the finger growth regime of convection mixing is particularly important. With the development of the visualization and surface measurement techniques, X-ray computed tomography (CT) [28] and Magnetic Resonance Image (MRI) are the available experimental methods to observe multiphase flow behavior in porous media [29]. MRI technology has the advantage of nondestructive visualizing the fluid distribution during fluid flow in porous media [30]. Up to date, few studies was conducted to experimental investigate the density-driven natural convection process in porous media. In this paper, we visualized the density-driven natural convection in porous media by using MRI. The development of convective fingers was investigated, and the quantitative analysis of the finger growth rates, the wave numbers, the initial wave length, the convection onset time and mixing time were achieved. And then we discussed the effect of Rayleigh number on the wave number, wave length and convection time. Such investigation aid in conceptual understanding of the density-driven natural convection and help inform studies for CO2 sequestration in saline aquifers. 2. Experimental methods 2.1. Experimental apparatus and material The schematic diagram of the experimental setup is shown in Fig. 1, including a high-resolution MRI machine for images acquisition, and a vessel for porous media. The Varian 400 MHz NMR system (Varian Inc., Palo Alto, CA, USA), 9.4 T magnet, was used for Proton (1H) images acquired. Pack sand was contained in a vessel with 150-mm length, 30mm internal diameter. The vessel was made from melamine resin which did not affect the magnetic signals during measurement and

was transparent to monitor outside MRI system as shown in Fig. 1. After degassing step in a Nalgene transparent polycarbonate desiccator (5311-0250, Thermo Scientific, Waltham, MA, USA) using a vacuum pump (Edwards, model E2M1.5), fluid was injected into the porous media by using an accurate syringe pump (Chin-Fin Technology, LC 3060). Four kinds of porous media were packed with plastic sands (US1, 250–425 μm in diameter, 337.5 μm in average diameter, Ube Sand Engineering Co. Ltd., XH40/60), and glass beads (BZ02, 150–250 μm, 200 μm in average diameter; BZ04, 350–450 μm, 400 μm in average diameter; BZ06, 550–650 μm, 600 μm in average diameter, As-One Co. Ltd., Japan). Three different kinds of fluid pairs: H2O/D2O, brine/D2O and H2O/ MnCl2 were used.D2O with a stated purity of 99.9% and the MnCl2 water solution were used as the denser fluids; deionized water or brine was used as the light fluid. The additive of coloured dye Bromocresol green (0.05%) was added into the lighter fluid to allow visual monitoring of flow and transport before MRI experiments. The properties of the porous media .i.e., porosity and permeability, and fluids were list in Table 1. Deionized water (H2O) and brine could be determined from the measured MRI intensities because 1H; However, deuterium oxide (D2O) does not emit an MR signal. The additive paramagnetic ion, such as Mn2 +, can reduce the relaxation time [31] and thereby the MnCl2 water solution is invisible in the MRI images. In other words, for the fluids system of H2O/D2O and brine/D2O, we can distinguish between 1H and D; for the fluid system of H2O/MnCl2, we can make a distinction between H2O and MnCl2 solution. 2.2. Experimental procedure After leakage check under the certain pressure, the vessel was filled with packed porous media. The light fluid (H2O or brine) was injected to

Table 1 Parameters of fluids and porous media in experiments. Case A B C D E F G H I J K L

Fluid pair H2O/D2O Brine/D2O H2O/MnCl2 H2O/D2O Brine/D2O H2O/MnCl2 H2O/D2O Brine/D2O H2O/MnCl2 H2O/D2O Brine/D2O H2O/MnCl2

Porous media

Porosity

Permeability

US1

0.354

28.63

BZ-02

0.48

57.1

BZ-04

0.398

64.64

BZ-06

0.431

111.58

Density difference (kg/m3) 0.1 0.05 0.1 0.1 0.05 0.1 0.1 0.05 0.1 0.1 0.05 0.1

Viscosity (Pa·s) 0.914 0.914 1.62 0.914 0.914 1.62 0.914 0.914 1.62 0.914 0.914 1.62

Diffusion coefficient (m2/s) −9

2 × 10 1.7 × 10−9 1.04 × 10−9 2 × 10−9 1.7 × 10−9 1.04 × 10−9 2 × 10−9 1.7 × 10−9 1.04 × 10−9 2 × 10−9 1.7 × 10−9 1.04 × 10−9

Ra 320.803 188.707 349.076 471.928 277.605 513.521 644.253 378.972 701.033 1027.05 604.150 1117.574

170

Y. Teng et al. / Magnetic Resonance Imaging 39 (2017) 168–174

Fig. 2. Longitudinal planes MRI images of the finger development of Brine/D2O fluid system, (a) case B of BZ02, (b) case H of BZ04, (c) case K of BZ06.

fully saturate porous media by using the vacuum chamber; and then denser fluid (D2O or MnCl2) was injected upward via the syringe pump at ambient conditions. The bottom connector tube was blocked carefully with the plug till the fluid reached the half height of the vessel. Before MRI experiments, we observed the convection process by visual inspection outside MRI first. Because the natural convection is caused immediately after the vessel turned upside down, we have to make sure the onset time (defined below) for finger forming is more than the MRI image acquisition time (2 min 8 s per image). Several experiments were repeated the above steps to guarantee there is enough for MRI image acquisition. After the calibration, one MRI scan was taken for the original state of fluids. Then take the vessel out of MRI, turn upside down quickly and re-place into the magnet. MRI images were simultaneously taken without interruption. MRI images were conducted using a spin echo multi-slice pulse sequence (SEMS), and the experimental parameters are as follows: TR = 1000 ms, TE = 5.62 ms, matrix size = 128 × 128, field of view (FOV) is 40 × 40 mm2 with 4.0 mm thickness, and the acquisition time is 2 min 8 s. As it is seen from Fig. 1, the interface of the two fluids was in the FOV. Images were taken in longitudinal planes and axial section. The scan was continued until variations in intensity of FOV was b1%. 2.3. Analytical parameters and dimensionless groups The experiments were conducted in four packed beds and three fluid systems to produce different Rayleigh number (Ra). Ra, as an dimensionless number, was used to indicate the onset of free convection and the degree of instability in porous media [32]. Ra is defined as k  Δρ  g  H Ra ¼ φμ D

where k is the permeability, Δρ is the density difference, g is acceleration due to gravity, H is the characteristics length of the system, ϕ is the porosity, μ is the viscosity and D is the diffusion coefficient [33]. Theoretical studies showed that the convection become preponderant when Ra N 4π2 [34,35]. All experiments in this study were designed to meet this criterion. Experiments specifications of fluid physical properties and Rayleigh number are reported in Table 1. The initial wavelength of the convective fingers were scaled with the inverse of Ra, and the critical time at which convective fingers occurs was inversely related to the square of Ra. the wavelength is given by [32] λ ¼ C1

μ φD Δρ  g  K

where the constant parameter C1 = 96.23 [32], λ is in meters. Linear stability analysis suggests the critical time should scale as [24]

T onset ¼ C 2

pffiffiffiffi!2 μ φ D Δρ  g  K

ð3Þ

C2 is the numerical constant in a wide range [24]; here, we choose C2 = 48.7 for the theoretical minimum of onset time. The factor C1 and C2 varies somewhat according to which method is used in the liner stability analysis. The time required to convective mixing exclude the diffusion dissolution time is calculated as [36] T mixing ¼

ð1Þ

ð2Þ

Hμ Δρ  g  K

Fig. 3. Axial section MRI images of the finger development of case K.

ð4Þ

Y. Teng et al. / Magnetic Resonance Imaging 39 (2017) 168–174

171

Fig. 4. Longitudinal planes MRI images of the one-finger pattern development, (a) case E of Brine/D2O in plastic sand, (b) case A of H2O/D2O in BZ02, (c) case C of H2O/MnCl2 in BZ02.

After the convective fingers occur, the time for fingers growth downwards is defined as mixing time Tmix. 3. Results and discussion Original longitudinal planes MRI images were cropped into 30 mm in width and 40 mm in height; and original axial section MRI images were cropped into 30 mm in width and 30 mm in height, those procedures were used to avoid redundant black background.

down, the interfaces changes to be uneven reflecting the occurrence of the fingers. The number and thickness of the fingers were different with porous media and fluids. At the time 62 min 2 s (case B), a one-finger pattern appeared on the interface and extend straight downward with time. When BZ04 glass beads were used as porous media (case H), the contour surfaces of the fingers became rough, and two convective fingers of different size could be recognized at 14 min 56 s in Fig. 2b. Four fingers with similar size occur on the interface in Fig. 2c can be seen from the interface at 6 min 24 s, Here we define the moments of finger appeared as the onset time of convection.

3.1. Convective finger growth regime and finger growth rates The longitudinal MRI images of the density-driven natural convection were shown in Fig. 2. The dark means the pore spaces filled with denser fluid (D2O or MnCl2 solution) and particles; in the contrast, bright means other pore spaces filled with H2O or brine. After the denser fluids injected in the vessel, the color was uniform and the interface between the denser (dark color) and lighter fluids (bright color) was flat, seen the left images as the original status. As the vessel turn upside

Fig. 5. Finger growth rate versus development path distance for one-finger pattern.

3.1.1. Convective finger growth regime The axial section images were simultaneously taken to analyze convective finger growth regime. Taking case K as an example, Fig. 3 shows the convective finger growth in slice 2. The convective fingers can be seen at 4 min 16 s in axial plane (Fig. 3 Frame 2), which was earlier than that from longitudinal planes. As time passed, the convective fingers penetrated downward with the number and thickness changing. The preferred finger growth was obvious and then propagated accompanying with mixing into larger finger. The red circles in Figs. 2c and 3 showed the coalescence of fingers due to small scale fingers that can merge adjacent fingers together to form fewer larger fingers before reaching the bottom. Fingers coalesce promote the wave and mixing which results in increase of wave length and mixing rate. Once the distance between the primary fingers was large enough, the finger coalescence stopped and the new fingers created. However, in several places, the new finger was usually swept by the light fluid that flow upward. In fact, in the beginning of mixing, small finger continuously appeared and grew as a bulge on the boundary layer between primary existing ones, seen in the red square in Fig. 2c. Overall, convective fingers growth regime may be divided into the following several stages approximately: convection onset period, the finger generation and propagate, the finger coalescence and new finger development. To investigate the flow dynamic properties simplicity, we focused on the properties of the single finger during convective development. As instances, the acquired MRI images of different convection process indicated that a single finger occurred in the cases A, B, C and E, as shown in Figs. 2a and 4. The finger development areas (the bottom half area of the FOV), was divided into 6 equal portions, as seen in Fig. 4b. According to the time that the fingertip reaches the certain position or the distance,

172

Y. Teng et al. / Magnetic Resonance Imaging 39 (2017) 168–174

Fig. 6. Axial section MRI images of at a particular time of Brine/D2O fluid system, (a) case B of BZ02 at 89 min 36 s, (b) case H of BZ04 at 21 min 20 s, (c) case H of BZ06 at 12 min 48 s.

we calculated the average growth rates (Fig. 5). The results revealed that the finger average growth rates do not remain the same at all position, reflecting three main regions with rising (A–C), stable (C–D), and declining (D–E) successively, and the slowest rate occurred at around the interface (A–B), compared with later regions. The reason was that diffusion dominates the convective at this region with the inconspicuous finger. At once finger appeared after onset time, the finger growth rate increased sharply at the region of B–C, because of the high density/concentration difference. However, with the density difference decreasing, the finger average growth rate tended to be constant (region of C–D). After that, a reduction in the growth rate occurred after the fingers reach the bottom region, because the denser fluid density was descended along the finger growth direction. 3.1.2. Effect of Rayleigh number on finger growth rate From Fig. 5, the finger growth rates varied with different Rayleigh number. It can be seen that during the process of finger extending downward, the rate increased with increasing of Ra, except case C. Ra was improved by increasing the permeability (case E) or density difference (case A) that compared with the minimum of case B; higher permeability or higher density difference both conducive to convection mixing. The density-driven convection only exist for case C; however, there was an isotopic exchange reaction H2O + D2O = 2HDO [37]

Fig. 7. Experimental measurements results and theoretical value of wave number (black) and initial wave length (red) versus Raleigh number. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

occurred in the fluid system containing D2O, which would enhance the fluid transports. This could explain why the finger growth rate of H2O/MnCl2 fluid system (case C) does not increase with the Rayleigh number. 3.2. Wave number and wave length Fig. 6 illustrated the finger distribution in various axial section at a particular time when the fingers had well-developed under different Ra conditions that corresponding to Fig. 2. The parts between the two lines in Fig. 6 specify the longitudinal planes location in axial section images. As a following procedure, the wave number is counted from the longitudinal planes images. Those parts with white color in Figs. 2 and 6 characterized the distribution of wave. The wave number of case B, case H and case K was distinct different in axial section MRI image. There's an obvious black round area in different slice along the finger growth direction for case B as shown in Fig. 6a. With the wave number rising, brine distributions in axial section got to be disconnected with an irregular shape. And also with increase of particle diameter, the convective fingers are more dispersedly distributed in the tube axial section (Fig. 6b and c). Also, we vary the Ra from 188 to 604 by using three kinds porous media (BZ02, BZ04 and BZ06) and fluids (density difference) system. The wave number were clearly influenced by the permeability of the

Fig. 8. Change in the normalized MRI intensity with time for H2O/MnCl2 fluid system.

Y. Teng et al. / Magnetic Resonance Imaging 39 (2017) 168–174

Fig. 9. Experimental measurements results of onset time by signal intensity and image sequence plotted against the theoretical value versus Rayleigh number.

porous media. With the Ra increase, it seems to be an upward trend for the wave number. The scaling relationship obtained from linear stability can be expressed as wave number a = 0.05 Ra. In order to verify the relationship between wave number and Ra, statistical analyses were applied to all cases to provide a quantitative interpretation. The experiments results of wave number as a function of Ra were compared with those obtained by linear stability analysis, as shown in Fig. 7. The scaling relationship obtained from MRI images (black dots) had the same changing tendency as the theoretical curve of linear stability analysis (black line). Results demonstrated that, the wave number increased with the increase Ra. The fingers have a diameter of approximately 5–25 mm (i.e. Figs. 2 and 4), which are orders of magnitude larger than particle diameters (0.2–0.6 mm). The linear stability analysis showed that the wave length of the fastest growing finger is defined as Eq. (2). The quantitative wave length was plotted as a function of Ra, as depicted in Fig. 7 in red color dots, the finger diameter was slightly smaller than the wave length estimated by linear stability analysis. All experiments data points followed a similar trend that the wave length was inversely proportional to the Ra. Two distinct stages can be identified in Fig. 7. For Ra b 500, the wave length of the initial convection instabilities varied as Ra− 1; while for Ra N 500, the normalized wave length was approximately constant λ ≈ 4 mm. In the last stage, the wave length showed no dependence on the density difference. The experimental results were consistent with those theoretical models supporting their validity. However, experimental value of wave length were less than the theoretical ones in Ra; for example, Ra = 188.707, the linear stability analysis predicted the wave length of 37.7 mm compared to an experimental value of 25 mm. This was due to the fact that the sample vessel limited the finger growth, the tube diameter was smaller than the one-finger pattern wave length. The design of the vessel was restricted by the MRI imaging probe diameter (40 mm).

173

Based on the MRI images, we estimated the normalized MRI intensity within the packed beds, as shown in Fig. 8. The MRI intensity was normalized until the signal tended to be stable for each experiment. The convection mixing process became fast with the increasing permeability of porous media. Using the signal intensity method, the onset time and the mixing time were measured as follows. The onset time was evaluated from the distance from the natural convection started (t = 0) to the intersection point of the regression line and normalized intensity equal 1. The mixing time was evaluated by subtracting the onset time from the time when the intersection point of the regression line and normalized MRI intensity reduces to 0 as shown in Fig. 8. Analytical and numerical calculations have been performed to predict the onset time and mixing time for density-driven natural convection by Eqs. (3) and (4), respectively. Fig. 9 shows the onset time achieved by two methods plotted against the theoretical value in Eq. (3). Fig. 10 shows the mixing time plotted against the theoretical value in Eq. (4). The experimental results of two methods were approximate; and both for onset time and mixing time were close to the theoretical values. The convection onset time corresponded to the end of diffusive mixing period when the fingers occurred, and the mixing time corresponded to the to the end of convection mixing period when the fingers became blurred. With the Ra increase, it seemed to be a downward trend for the convective fingers occurs time. In other words, the convective fingers occurred quickly with an increasing in Ra. The variation of mixing time to Ra was quite similar with the onset time, the mixing time was reduced with raising of Ra. The field of view (FOV) of our experiment was only 40 mm in height, but the length of the vessel for convection mixing is 150 mm. When the theoretical value of mixing time was calculated to comparing with the experimental results, H is the length of FOV instead the length of the vessel. Because the limitation of FOV, the experimental results of mixing time was shorter than onset time. In some ways, the onset and mixing time showed the same changing rule as the wave length; for Ra N 500, the time value reached a plateau, onset time and mixing time not change with the Ra. All those results confirm the dominant role of Ra as the governing dimensionless parameter in density-driven natural convection process. However, there are also fluctuations in the relation curve between time and Ra. Those fluctuated position in the curves was outlined in a dotted box as shown in Figs. 9 and 10. The fluid system of those cases with fluctuations all is H2O/MnCl2, a possible explanation could be that, as we have introduced in the before section, there was only density-driven convection occurred in H2O/MnCl2 fluid system; compared with other fluid systems (i.e. H2O/D2O), the H2O/MnCl2 fluid system is lack of isotopic exchange reaction; therefore, no matter

3.3. Critical onset time and mixing time of miscible fluids In order to evaluate the time of the convective mixing process, we measured the critical time and mixing time through two methods: image sequence and signal intensity. The MRI image acquisition time was 2 min 8 s and image capturing was uninterrupted, it was convenient to determine the onset time and mixing time according to the MRI images. As shown in Figs. 2–4, every image below was labeled the end of acquisition time. The disadvantage to this approach was the accurate for time detecting, the point in time was always a multiple of 2 min 8 s (128 s).

Fig. 10. Experimental measurements results of mixing time by signal intensity and image sequence plotted against the theoretical value versus Rayleigh number.

174

Y. Teng et al. / Magnetic Resonance Imaging 39 (2017) 168–174

diffusive mixing period or convection mixing period, the H2O/MnCl2 fluid system take more longer time. 4. Conclusion The setup in this study provided an analog for CO2/brine densitydriven natural convection in a vertical tube packed with porous media. The study focused on density-driven flow transport processes that occur in a vertical tube. A high-resolution MRI system captured the evolution of convective fingers with time from which key regimes of the fingers growth were determined. Four main regimes were found in the convective finger growth process: convection onset period, the finger generation and propagate, the finger coalescence and new finger development. While for a single finger case, the finger average growth rates changed with time and position, reflecting three main regions with rising, stable, and declining successively, and the slowest rate occurred at around the interface because of the diffusion dominate. Also the finger growth rate and the wave number increased with an increase of Rayleigh number; however, the wave length was inversed proportion to Rayleigh number and smaller than the theoretical value due to the limited of tube diameter. The onset time and mixing time of convective fingers decreased with the Rayleigh number increase. The wave length, onset time and mixing time kept to be constant until Rayleigh numbers were larger than 500. The work provided a new idea and method for reducing dissolution trapping time of CO2 sequestration. Acknowledgment This study has been supported by the National Natural Science Foundation of China (Grant No. 51506024), the National Key Research and Development Program of China (Grant No. 2016YFB0600804). It has been also supported by the Fundamental Research Funds for the Central Universities (DUT16QY12). Authors appreciated the help of Prof. Suekane from Tokyo Institute of Technology. References [1] Khudaida KJ, Das DB. A numerical study of capillary pressure–saturation relationship for supercritical carbon dioxide (CO2) injection in deep saline aquifer. Chem Eng Res Des 2014;92(12):3017–30. [2] Nazari Moghaddam R, et al. Quantification of density-driven natural convection for dissolution mechanism in CO2 sequestration. Transp Porous Media 2011;92(2): 439–56. [3] Kneafsey TJ, Pruess K. Laboratory flow experiments for visualizing carbon dioxideinduced, density-driven brine convection. Transp Porous Media 2009;82(1):123–39. [4] Pau GSH, et al. High-resolution simulation and characterization of density-driven flow in CO2 storage in saline aquifers. Adv Water Resour 2010;33(4):443–55. [5] Khosrokhavar R, et al. Visualization and investigation of natural convection flow of CO2 in aqueous and oleic systems. J Pet Sci Eng 2014;122:230–9. [6] Riaz A, et al. Onset of convection in a gravitationally unstable diffusive boundary layer in porous media. J Fluid Mech 2006;548:87–111. [7] Javaheri M, Abedi J, Hassanzadeh H. Linear stability analysis of double-diffusive convection in porous media, with application to geological storage of CO2. Transp Porous Media 2010;84(2):441–56. [8] Backhaus S, Turitsyn K, Ecke R. Convective instability and mass transport of diffusion layers in a Hele-Shaw geometry. Phys Rev Lett 2011;106(10):104501.

[9] Tsai PA, Riesing K, Stone HA. Density-driven convection enhanced by an inclined boundary: implications for geological CO2 storage. Phys Rev E 2013;87(1):011003. [10] Weir GJ, White SP, Kissling WM. Reservoir storage and containment of greenhouse gases. Transp Porous Media 1996;23(1):37–60. [11] Lindeberg E, Wessel-Berg D. Vertical convection in an aquifer column under a gas cap of CO2. Energy Convers Manag 1997;38:S229–34. [12] Riaz A, et al. Onset of convection in a gravitationally unstable diffusive boundary layer in porous media. J Fluid Mech 2006;548(1):87. [13] Cooper CA, Glass RJ, Tyler SW. Experimental investigation of the stability boundary for double-diffusive finger convection in a Hele-Shaw cell. Water Resour Res 1997; 33(4):517–26. [14] Kneafsey TJ, Pruess K. Laboratory flow experiments for visualizing carbon dioxideinduced, density-driven brine convection. Transp Porous Media 2010;82(1):123–39. [15] Mojtaba S, et al. Experimental study of density-driven convection effects on CO2 dissolution rate in formation water for geological storage. J Nat Gas Sci Eng 2014;21: 600–7. [16] Thomas C, et al. Experimental study of CO2 convective dissolution: the effect of color indicators. Int J Greenhouse Gas Control 2015;42:525–33. [17] Wylock C, et al. Experimental study of gas–liquid mass transfer coupled with chemical reactions by digital holographic interferometry. Chem Eng Sci 2011;66(14): 3400–12. [18] Faisal TF, Chevalier S, Sassi M. Experimental and numerical studies of density driven natural convection in saturated porous media with application to CO2 geological storage. Energy Procedia 2013;37:5323–30. [19] Neufeld JA, et al. Convective dissolution of carbon dioxide in saline aquifers. Geophys Res Lett 2010;37(22) n/a-n/a. [20] Yang C, Gu Y. Accelerated mass transfer of CO2 in reservoir brine due to density-driven natural convection at high pressures and elevated temperatures. Ind Eng Chem Res 2006;45(8):2430–6. [21] Farajzadeh R, et al. Mass transfer of CO2 into water and surfactant solutions. Pet Sci Technol 2007;25(12):1493–511. [22] Kneafsey TJ, Pruess K. Laboratory experiments and numerical simulation studies of convectively enhanced carbon dioxide dissolution. Energy Procedia 2011;4: 5114–21. [23] Walker KL, Homsy GM. Convection in a porous cavity. J Fluid Mech 1978;87(03): 449–74. [24] Lindeberg E, Wessel-Berg D. Upscaling studies of diffusion induced convection in homogeneous and heterogeneous aquifers. Energy Procedia 2011;4:3927–34. [25] Fernandez J, et al. Density-driven unstable flows of miscible fluids in a Hele-Shaw cell. J Fluid Mech 2002;451:239–60. [26] Hidalgo JJ, et al. Scaling of convective mixing in porous media. Phys Rev Lett 2012; 109(26):264503. [27] Emami-Meybodi H, Hassanzadeh H. Two-phase convective mixing under a buoyant plume of CO2 in deep saline aquifers. Adv Water Resour 2015;76:55–71. [28] Wang L, et al. Three-dimensional visualization of natural convection in porous media. Energy Procedia 2016;86:460–8. [29] Song Y, et al. Magnetic resonance imaging study on near miscible supercritical CO2 flooding in porous media. Physics of Fluids (1994-present) 2013;25(5):053301. [30] Song Y, et al. An experimental study on CO2/water displacement in porous media using high-resolution Magnetic Resonance Imaging. Int J Greenhouse Gas Control 2012;10:501–9. [31] Bloembergen N, Purcell EM, Pound RV. Relaxation effects in nuclear magnetic resonance absorption. Phys Rev 1948;73(7):679–712. [32] Xu X, Chen S, Zhang D. Convective stability analysis of the long-term storage of carbon dioxide in deep saline aquifers. Adv Water Resour 2006;29(3):397–407. [33] Heidorn S-C, et al. Consecutive mechanism in the diffusion of D2O on a NaCl (100) bilayer. ACS Nano 2015;9(4):3572–8. [34] Lapwood E. Convection of a fluid in a porous medium. Mathematical Proceedings of the Cambridge Philosophical Society. Cambridge Univ Press; 1948. [35] Combarnous M. Hydrothermal convection in saturated porous media. Adv Hydrosci 1975;10:231–307. [36] Suekane T, Sakai S. Three-dimensional visualization of natural convection of miscible fluids due to the density difference in a packed bed. 2015;1:475–86. [37] Pyper JW. Equilibrium constants of hydrogen–deuterium–tritium self-exchange reactions in water vapor as studied with a pulsed molecular-beam quadrupole mass filter. J Chem Phys 1975;62(7):2596.