Subscriber access provided by UNIV OF ARIZONA
Article
Dissolution Trapping of Carbon Dioxide in Heterogeneous Aquifers Mohamad Reza Soltanian, Mohammad Amin Amooie, Naum Gershenzon, Zhenxue Dai, Robert Ritzi, Fengyang Xiong, David R. Cole, and Joachim Moortgat Environ. Sci. Technol., Just Accepted Manuscript • Publication Date (Web): 09 Jun 2017 Downloaded from http://pubs.acs.org on June 12, 2017
Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.
Environmental Science & Technology is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.
Page 1 of 33
Environmental Science & Technology
Dissolution Trapping of Carbon Dioxide in Heterogeneous Aquifers Mohamad Reza Soltanian,†,‡ Mohammad Amin Amooie,† Naum Gershenzon,¶ Zhenxue Dai,§,k Robert Ritzi,¶ Fengyang Xiong,† David Cole,† and Joachim Moortgat∗,†,⊥ †School of Earth Sciences, The Ohio State University, Columbus, OH 1
‡
[email protected], Exponent, 1055 East Colorado Boulevard, Suite 500 Pasadena, CA 91106 ¶Department of Earth and Environmental Sciences, Wright State University, Dayton, OH §Los Alamos National Laboratory, Los Alamos, New Mexico, USA kCollege of Construction Engineering, Jilin University, Changchun, China ⊥Correspondence: 303 Mendenhall Laboratory, 125 South Oval Mall, Columbus, OH 43210, Tel: +1-614-688-2140, Fax: +1-614-292-7688 E-mail:
[email protected] 2
Abstract
3
The geologic architecture in sedimentary reservoirs affects the behavior of density-
4
driven flow and the dispersion of CO2 -rich brine. The spatial organization and con-
5
nectivity of facies types play an important role. Low-permeability facies may sup-
6
press fingering and reduce vertical spreading, but may also increase transverse mix-
7
ing. This is more pronounced when geologic structures create preferential flow path-
8
ways through connected facies types. We perform high-resolution simulations of three-
9
dimensional (3D) heterogeneous formations whose connectivity cannot be represented in
1
ACS Paragon Plus Environment
Environmental Science & Technology
10
two-dimensional models consistent with percolation theory. This work focuses on the
11
importance of 3D facies-based heterogeneity and connectivity on advection-diffusion
12
transport of dissolved CO2 . Because the dissolution of CO2 and the subsequent density
13
increase of brine are the driving force for gravitational instabilities, we model the phase
14
behavior with the accurate cubic-plus-association equation-of-state, which accounts for
15
the self-association of polar water molecules and the cross-association between CO2 and
16
water. Our results elucidate how the spatial organization of facies affects the dynamics
17
of CO2 convective mixing. Scaling relations for the evolution of a global dispersion-
18
width provide insights that can be universally applied. The results suggest that the
19
long-term evolution and scaling of dispersion are surprisingly similar for homogeneous
20
and (binary and multiscale) heterogeneous porous media.
21
Introduction
22
Carbon dioxide (CO2 ) injection in geologic formations is a promising short- and long-term
23
strategy to slow the increase in atmospheric CO2 concentration. 1–5 To assure decision makers
24
and the public of the feasibility, economics, and storage permanence of Carbon Capture
25
and Storage (CCS), a thorough understanding of the basic science behind CCS is critical.
26
However, comprehensive assessment of CCS is challenging due to the complex physical and
27
chemical processes controlling the transport of CO2 in geologic reservoirs. This includes
28
CO2 phase behavior and different types of trapping mechanisms (e.g., residual, dissolution,
29
structural, and mineral). 6 In this work, we focus on dissolution or solubility trapping of CO2
30
in brine, which is believed to be an intermediate-term (0–1000 years) mechanism. The focus
31
is on the dissolution mechanisms taking place in the aquifer after a plume reaches the cap
32
rock. Residual, 7 structural, 8 and mineral trapping are not considered.
33
Solubility trapping is enhanced when gravitational instabilities (fingering) are triggered
34
by a local increase in brine density as CO2 dissolves into brine in the top of an aquifer. 9–13
35
Whether the interface between CO2 -rich brine and fresh brine becomes gravitationally un2
ACS Paragon Plus Environment
Page 2 of 33
Page 3 of 33
Environmental Science & Technology
36
stable depends on the ratio of advection to molecular diffusion. When a fingering instability
37
is triggered, dissolved CO2 is mixed throughout the aquifer at advective time-scales, which
38
can be much faster than diffusive transport. This fast spreading towards higher depths im-
39
proves the storage capacity of a given aquifer and decreases the leakage risk in case of cap
40
rock failure.
41
Heterogeneity in petrophysical properties (e.g., permeability and porosity) has a signif-
42
icant impact on density-driven flow of CO2 . 11,13–15 Both rock and fluid properties directly
43
affect the characteristics of fingers including the critical time for the onset of the instability
44
(tc ) and a critical wavelength of fingers (λc ). The results of linear stability analyses, as well
45
as numerical simulations show that: 16–21
46
47
φµ , k∆ρg 2 φµ = c2 D , k∆ρg
λc = c 1 D
(1)
tc
(2)
48
where µ is brine viscosity, φ is porosity, D is the diffusion coefficient for CO2 in brine, k is
49
permeability, ∆ρ is the maximum density increase of the aqueous phase upon CO2 disso-
50
lution, and g is the gravitational acceleration. The c1 and c2 are proportionality constants
51
that vary by orders of magnitudes in the cited references.
52
Convective mixing and solubility trapping can be enhanced by subsurface heterogene-
53
ity. 11,14 Geological heterogeneity is controlled by the spatial organization of facies types which
54
are three-dimensional (3D) bodies of rocks, whose differentiation provides a useful concep-
55
tual framework to characterize heterogeneity in attributes of interest (e.g., permeability). 22
56
The architecture of facies types may exhibit sharp discontinuities across boundaries between
57
depositional features, for instance across the sandstone-shale contacts that are common in
58
fluvial deposits. In fact, some potential CO2 reservoirs have sedimentary architecture re-
59
flecting fluvial deposition with complex architecture created by depositional processes. Such
60
processes often result in variable connectivity and highly tortuous flow pathways. 6,23–25 Well-
3
ACS Paragon Plus Environment
Environmental Science & Technology
61
known examples include the Victor interval of the Ivishak Formation, 8,26 lower Mt. Simon
62
Sandstone, 25 the lower Tuscaloosa Formation, 23,27 the lower Paaratte Formation, 27,28 and
63
the Morrow Sandstone. 29
64
Facies connectivity significantly affects CO2 migration by creating preferential flow path-
65
ways. 8,13,30–32 Connected zones of high-permeability facies exhibit shorter critical times and
66
wavelengths, and thus experience more pronounced fingering. In this work we perform 3D
67
simulations to (1) represent 3D geologic structure and connectivity, (2) understand the fluid
68
dynamics of gravitational fingering when 3D heterogeneity is considered, and (3) provide
69
scaling relations to predict the long-term dissolution trapping of CO2 . Earlier 3D simulation
70
studies of CO2 convective mixing have only considered homogeneous media. 19
71
The volume proportion of each facies is the most important factor controlling the con-
72
nectivity. 33 Percolation theory suggests that randomly placed cells of one type have at least
73
one connected tortuous pathway through a cluster that spans opposing boundaries if the
74
volume proportion of that type is above the 3D percolation threshold of 31.16%. 33,34 If the
75
spatial organization of a facies type is represented within a finite domain, a tortuous fully
76
connected pathway is formed at even lower volume proportions. 35 In this work we quantify
77
connectivity by studying connected clusters of high-permeability cells within two different
78
types of geologic models using information on how geologic structure affects the percolation
79
principle. 35–38 In the first model we create synthetic binary distributions of facies types using
80
a Markov Chain approach. Facies architectures are represented by physically quantifiable
81
attributes such as mean lengths and volume proportions. The second geologic model is a
82
multiscale complex representation of a fluvial architecture, which consists of several facies
83
types spanning from the cm to the hundred meter scale. In both models, we create different
84
realizations, each with a different volume proportion of high-permeability facies.
85
Our main objective is to study the impact of three-dimensional facies connectivity on
86
the density-driven convective mixing of dissolved CO2 . By using models for which the con-
87
nectivity of facies types are known, we aim to quantify the role of connectivity on solubility
4
ACS Paragon Plus Environment
Page 4 of 33
Page 5 of 33
Environmental Science & Technology
88
trapping in general, and the dynamics of CO2 mixing in particular.
89
Methodology
90
Heterogeneity structures and connectivity
91
We consider (1) a synthetic binary heterogeneity model, and (2) a heterogeneity model that
92
reflects a multiscale fluvial architecture. In both models we vary the volume proportion
93
of high-permeability facies to study connectivity. There are several approaches to study
94
connectivity including the use of (1) flow- and transport-based measures, 39 (2) critical path
95
analysis, 40 and (3) facies-based approaches linked with percolation theory. 8,35,36 Here we use
96
the facies-based approach and define the cell connectivity in three dimensions.
97
Binary heterogeneity
98
We use generic heterogeneity models in binary media with spatial distribution of high-
99
permeability facies within a background of lower permeability facies. We model a 50 m ×
100
50 m × 50 m domain, discretized by a fine 80 × 80 × 80 grid with 512,000 grid cells.
101
We show the fundamental importance of 3D heterogeneities arising from facies archi-
102
tecture. We choose a parsimonious set of two facies types. 41,42 The multifacies model is
103
described in the next section. Different structures are created using a Markov Chain model
104
and indicator simulation with quenching in T–PROGS. 43 Rather than focusing on any sin-
105
gle site or data set, we represent general permeability distributions consisting of high- and
106
low-permeability facies. The volume proportions of high-permeability facies are increased
107
from 5 to 16, 28, and 40% (Figure S1 in the Supporting Information shows realizations with
108
16 and 28%). In all binary cases, we use an isotropic integral scale (λ) of 10 m. The effect of
109
small-scale anisotropy within grid cells is presented as a special case in the Results and Dis-
110
cussion. A high and a low permeability of 5, 000 md and 500 md are assigned to the two facies
111
types, respectively. The porosity of high- and low-permeability facies are calculated as 0.34 5
ACS Paragon Plus Environment
Environmental Science & Technology
112
and 0.13, respectively, using k = aφb where a and b are constants. 44,45 We do not consider
113
permeability variations within each facies but assume that the spatial organization of facies
114
types is the major source of heterogeneity. We choose domain size, cell size, and permeability
115
values based on prior experience and Equation (1), which shows that λc is proportional to
116
1/k. Thus, the cell size should be small enough to capture the smallest fingers within high
117
permeability facies. At the same time, the domain size should be large enough compared to
118
λc so that the domain encompasses the largest fingers within low-permeability facies. For
119
comparison, we also create homogeneous domains with permeabilities and porosities based
120
on the means of the heterogeneous domains (Table S1).
121
We use a search algorithm to determine how high-permeability facies are connected in
122
the 3D heterogeneous binary model. As an example, for the volume proportion of 28% this
123
algorithm determined that 96% of high-permeability facies (27% of all cells) are connected
124
in a single cluster that spans between any pair of opposing domain boundaries, creating
125
a percolating cluster (see Figure S1d). It is obvious that the volume proportion of 28%
126
is below the volume proportion of 31.16% known as percolation threshold for a random,
127
infinite cubic lattice. However, the structure added to the finite-domain model reduced the
128
threshold, as is consistent with previous studies. 35,38,46 The percentage of connected cells
129
for volume proportions of 5, 16, and 40% are ≪ 1, 1.4, and 39.6%, respectively. There
130
is no connected percolating cluster for volume proportions of 5 and 16%. However, there
131
is a connected percolating cluster for volume proportions of both 28 and 40%. Note that
132
two-dimensional (2D) slices from a 3D model do not have the connectivity of the 3D domain
133
they are extracted from, as per percolation theory. 46
134
Multiscale heterogeneity
135
Many candidate reservoirs for geologic CO2 sequestration have fluvial depositional structures.
136
Recent studies have led to new conceptual and quantitative models for such architectures
137
over a range of spatial scales that are relevant to fluid flow in oil and deep brine reservoirs. 47
6
ACS Paragon Plus Environment
Page 6 of 33
Page 7 of 33
Environmental Science & Technology
138
This approach to characterizing fluvial architectures, including scales of stratification, volume
139
proportions, relative sizes, and 3D geometry for facies types at each scale, was incorporated
140
into a multiscale heterogeneity model by Guin et al. 38 and Ramanathan et al. 37
141
The multiscale heterogeneity model 8,25,41,48 represents fluvial deposits dominated by sandy
142
gravel and the associated conglomeratic reservoirs. At the largest scale, the model comprise
143
active channels and compound bars. The compound bars, in turn, comprise unit bars.
144
Within the unit bars there are the smallest scale cross-strata with different textures. Im-
145
portantly, there are high-permeability sets of cross-stratified open-framework conglomerate
146
(OFC), which are known to create preferential flow pathways. The individual OFC cross-sets
147
are represented in the model, and clusters of continuously connected OFC cells create pref-
148
erential flow pathways. Guin et al. 38 found that tortuous spanning pathways exist through
149
connected OFC cells if OFC volume proportions were at or above 20%, independent of other
150
input parameters. The number, size, and orientation of such OFC clusters change with the
151
volume proportions of each facies. Even for a given volume proportion, these properties
152
change across hierarchical levels. 37,38 For instance, connected OFC cells at the smallest scale
153
(level I) form flow pathways that vertically span (i.e., connect across) a single unit bar at the
154
next scale (level II). Connections across unit bars at level II enhance lateral branching and
155
the many clusters within unit bars connect into a smaller number of larger clusters at the
156
scale of multiple unit bars. At the scale of a compound bar deposit (level III), these clusters
157
are typically connected into one or two large clusters that span all the domain boundaries.
158
We use this complex heterogeneity model to investigate the effect of OFC connectivity on
159
the dynamics of density-driven flow of CO2 .
160
We generate three realizations in which the volume proportions of OFC are 19, 24, and
161
28%. The volume proportion controls the percentage of OFC cells that are connected in
162
clusters to create preferential flow pathways. 8 Such spanning clusters may develop at OFC
163
proportions below the well-known percolation threshold of 31.16%. The search algorithm
164
described in the previous section is used to identify connected clusters of OFC cells, and
7
ACS Paragon Plus Environment
Environmental Science & Technology
165
to quantify their attributes. The proportion of all connected cells are 2.5, 17.4, and 25.6%,
166
respectively, for volume proportions of 19, 24, and 28% OFC cells. There are spanning
167
clusters in the cases with 24 and 28% OFC, but not for 19%. Therefore, the multiscale
168
heterogeneity models range from one with no spanning cluster (19%) to one with about 90%
169
connected OFC in one spanning cluster (28%).
170
The permeability fields of cases with 19 and 28% volume proportions are illustrated in
171
Figure S2. The domain size is 100 m × 100 m × 15 m, discretized by a 2 m × 2 m × 0.05 m
172
grid-cell size with 750,000 grid cells. For each facies type we use lognormal permeability
173
distributions as found in natural deposits. The permeability of each cell is assigned from
174
the probability density functions (PDFs) defined for each facies type (see Figure 3, panels
175
on the right, in Ramanathan et al. 37 ). The geometric mean of permeabilities for 28, 24, and
176
19% cases are 720, 891, and 1051 md, respectively.
177
The thickness-to-length ratio in the geometry of facies is significantly less than unity at
178
all spatial scales. Therefore, the heterogeneity gives rise to an effective anisotropy at the
179
scale of the sedimentary facies types. Anisotropy within the 5 cm thickness of grid cells
180
is not considered (as in Gershenzon et al. 8 ). Note that the grid-cell size is defined by the
181
characteristic size of the smallest scale facies types in the model architecture. Larger grid-
182
cell sizes cannot capture the small-scale heterogeneities. Therefore, to make the problem
183
computationally feasible, there is a trade-off between grid-cell size and domain size.
184
Numerical methods and boundary conditions
185
We use numerical methods presented in earlier works. 21,49–56 The governing equations are
186
presented in Supporting Information, Text S1. We consider a fluid mixture with two compo-
187
nents (water and CO2 ). Note that we consider fluid compression due to the pressure build-up
188
during CO2 injection into an aquifer with impermeable boundaries and volume swelling of
189
brine due to CO2 dissolution. This is in contrast to previous studies on gravitational fingering
190
in which an incompressible aqueous phase is considered. 8
ACS Paragon Plus Environment
Page 8 of 33
Page 9 of 33
Environmental Science & Technology
191
To solve the governing equations we use higher-order finite element methods. Specifi-
192
cally, we adopt a second-order discontinuous Galerkin (DG) method to solve the transport
193
equation. 57 Using DG is beneficial in studying flow instabilities since it is able to capture
194
the small-scale onset of viscous and gravitational fingers on relatively coarse grid cells by
195
reducing numerical dispersion. This helps us to resolve λc with only few grid cells. 21 We
196
use the mixed hybrid finite element (MHFE) method to simultaneously solve for continuous
197
pressure and velocity fields. MHFE provides an accurate and globally continuous velocity
198
field for highly heterogeneous domains. This feature is critical given the sharp contrasts in
199
the permeability fields arising from the spatial organization of sedimentary facies.
200
Phase behavior is based on the cubic-plus-association (CPA) equation-of-state (EOS).
201
The CPA-EOS was developed to model the phase behavior of fluids that contain polar
202
molecules, such as water 58 or asphaltenes, 59 by using thermodynamic perturbation theory.
203
The CPA-EOS takes into account the self-association between water molecules, as well as
204
the (polarity-induced) cross-association between water and CO2 . 51,60,61 This approach is an
205
improvement over previous studies that relied on empirical correlations for the brine density
206
and Henry’s law for CO2 solubilities. We strictly enforce local thermodynamic equilibrium.
207
Also, the CPA-EOS has been shown to accurately reproduce measured densities of CO2 -
208
brine mixtures and the volume swelling due to CO2 dissolution. Volume swelling and the
209
associated movement of the interface between free CO2 in a gas cap and the underlying brine,
210
cannot be modeled by considering only the aqueous phase 17 and replacing the gas cap with
211
a constant-CO2 -composition Dirichlet boundary condition. 62
212
The temperature and the pressure at the bottom boundary are, respectively, 77 ◦ C and
213
100 bar. At these conditions, the aqueous phase density is ρw = 981 kg/m3 , the maximum
214
CO2 solubility is 1.7 mol%, and the CO2 -saturated aqueous phase density is ρ = 987 kg/m3
215
with ∆ρ = (ρ − ρw )/ρw = 0.6% (note that CO2 behaves as a supercritical fluid in this
216
pressure and temperature condition). The diffusion coefficient is D = 1.33 × 10−8 m2 s−1 .
217
We consider 3D bounded domains with impermeable no-flow Neumann conditions for all
9
ACS Paragon Plus Environment
Environmental Science & Technology
Page 10 of 33
218
boundaries. The domain is initially saturated with brine. CO2 is injected uniformly into
219
the top grid cells at a low rate (as a source term of 0.1% pore volume (PV) per year). The
220
domain remains in a single aqueous phase without a capillary transition zone. 63,64 Because
221
of the impermeable boundaries the pressure increases upon injection. We also simulated
222
a case with an open bottom boundary in which the pressure remains roughly constant. 64
223
However, the evolution of the dispersion width (and the mixing dynamics more generally) do
224
not appear to be sensitive to this boundary condition. The reason is that, while the pressure
225
evolution is quite different for the two different boundary conditions, the density difference
226
within the aqueous phase upon CO2 dissolution is nearly identical, and this is the driver for
227
convective mixing. Additional testing shows that the computational domain is sufficiently
228
wide that the lateral boundaries do not impact the temporal evolution of CO2 mixing.
229
The density-driven flow of CO2 is quantified by the global dispersion-width σz (the square
230
root of the spatial variance σz2 ) of the molar density of CO2 following Soltanian et al. 13 (see
231
also Agartan et al. 14 ). The choice of CO2 molar density, rather than composition, accounts
232
for the density change due to both the CO2 dissolution and the compressibility of aqueous
233
phase. The first three central moments of the distribution define the spatial variance σz2 of
234
the molar density of CO2 (C ≡ czCO2 ): σz2 (t)
N002 = − N000
N001 N000
2
(3)
235
where N001 /N000 expresses the coordinate location of the center of mole in the vertical (z)
236
direction. The Nijk is the ijk-th moment of the CO2 molar density distribution in space,
237
analogous to the moment of concentration distribution, 65 defined here as:
Nijk (t) =
Z
Lx 0
Z
Ly 0
Z
Lz
xi y j z k φC(x, y, z, t)dxdydz
(4)
0
238
with Lx , Ly , and Lz the domain size in x, y, and z directions, respectively. The total number
239
of moles of CO2 in the domain is physically represented by the zeroth spatial moment, N000 . 10
ACS Paragon Plus Environment
Page 11 of 33
Environmental Science & Technology
240
The N000 does not remain constant over time since CO2 is continuously injected into the
241
domain. The spatial variance is not sensitive to N000 increase since the first and second
242
moments are normalized by N000 . The variances in the x-direction and y-direction are
243
negligible due to the dominant downward migration of CO2 -rich brine. Fundamentally σz is
244
a measure of spreading and dispersion in porous media, however, it may also represent the
245
mixing state since more spreading usually leads to more mixing. 66
246
Results and Discussion
247
Binary media
248
We first study the dynamics of fingering in the synthetic binary media and their homogeneous
249
equivalents.
250
We provide Figure 1a that visually shows the profiles for mole fraction of CO2 after 3, 4.5,
251
6.5, 8.5, 11, and 13.5 years for homogeneous domains corresponding to binary media with
252
16 and 28% volume proportions of high-permeability facie types. The mean (homogeneous)
253
permeability and porosity increase for larger volume proportions of high-permeability facies
254
(see Table S1). It is clear that the onset of fingering occurs later for the homogeneous
255
case of 16% than that of 28% and there are more CO2 fingers initiated for the 28% volume
256
proportion. Using the values in Table S1 we find that λc (28%)/λc (16%) is about 0.8 and
257
tc (28%)/tc (16%) is about 0.65.
258
Phase behavior is a critical part of our simulations. A thermodynamically consistent
259
description of phase behavior results in a driving force for gravitational instabilities (∝ ∆ρ)
260
that evolves non-trivially upon CO2 dissolution (see Figure 5 in Soltanian et al. 13 ). This
261
evolution of ∆ρ results in four distinct flow regimes in our scaling analyses of homogeneous
262
simulations presented in Figure 2a, which is consistent with our novel findings for 2D simu-
263
lations in homogeneous media. 13 The alternating flow regimes are interpreted as follows.
264
The first flow regime corresponds to the well-known diffusive transport of CO2 before the 11
ACS Paragon Plus Environment
Environmental Science & Technology
ACS Paragon Plus Environment
Page 12 of 33
Page 13 of 33
Environmental Science & Technology
265
onset of gravitational fingers (see Figure 1a, panel 3 years for 16%). As shown in Figure 2a,
266
the σz exhibits a classical Fickian scaling of σz ∝ t0.5 . The front in this regime travels to a
267
diffusion penetration depth, which obeys a (Dt)0.5 Einstein-type scaling.
268
The second flow regime corresponds to the onset of fingering and can be recognized as
269
the transition to a convection-dominated regime, according to the sharp increase in σz in
270
Figure 2a. The dependence of tc and λc on k and φ imply that higher permeabilities result in
271
an earlier onset time with a larger number of small-scale fingers, as illustrated in Figure 1a.
272
The initial scaling of σz in Figure 2a for the fingering regime is σz ∝ t3.5 . This is also
273
consistent with our 2D simulation results. 13 Note that the onset time is earlier for larger
274
values of homogeneous permeability and porosity.
275
The second regime is associated with reduction in ∆ρ within the fast sinking fingers. 13
276
This is because the high rate of downward advection leads to stretching of CO2 -rich plume,
277
hence lowering the peak CO2 concentration and subsequently the maximum ρ. Meanwhile,
278
the stretching-induced sharp concentration gradients across the finger boundaries cause the
279
diffusion of CO2 into the ambient upwelling brine. This results in depletion of CO2 in the
280
sinking fingers, which further contributes to ∆ρ reduction. In a negative feedback loop, the
281
resulting reduced ∆ρ inhibits further stretching of fingers (note that ∆ρ is the driving force).
282
The second regime is therefore transient, and is followed by a third restoration regime.
283
The third regime is attributed to a fingering stagnation period, in which the depleted
284
fingers become replenished with new CO2 that is injected from the top boundary and dis-
285
solved in brine at a diffusive rate. We find diffusive scaling in this regime, mainly because of
286
persistent coalescence and merging of growing fingers resulting in coarsened fingers, as well
287
as the transverse diffusion from the CO2 -rich fingers into the brine as the CO2 concentra-
288
tions inside the fingers become higher again during the restoration. This regime transitions
289
into the final fully-developed convection once the ∆ρ, through constant CO2 injection, again
290
exceeds a critical limit that is required to drive the gravitational instabilities. 13
291
Finally, in the fourth regime the flow is again advection-dominated and shows a sharp
13
ACS Paragon Plus Environment
Environmental Science & Technology
292
increase in the growth rates of σz . This regime and its related scaling for σz are the most
293
important ones since they control the long-term evolution of σz , while the first three regimes
294
occur at shorter time-scales (∼10 years in our results). What defines the long-term fate of the
295
injected CO2 is the evolution of this fourth flow regime, the extent of which is determined by
296
the vertical size of the aquifer. Thick aquifers may be more effective in trapping CO2 through
297
solubility. It has been also shown that the capillary sealing efficiency of cap rocks –another
298
storage mechanism– is enhanced by CO2 sequestration at higher depths. 67 This suggests
299
that CCS projects may benefit from deeper formations, in spite of the higher costs involved.
300
Off-shore sedimentary aquifers, which can have thicknesses of the order of a kilometer, 71 are
301
attractive for the same reasons.
302
Note that such distinct flow regimes, especially the second and third regimes, are not
303
observed in models that assume incompressible flow and a constant pressure and CO2 -
304
composition at the top boundary of a brine-saturated domain. Under such conditions, the
305
maximum value for the density contrast (∆ρmax ) between pure and CO2 -saturated brine
306
remains constant. A constant value of ∆ρmax is also assumed in all linear stability analyses
307
that lead to expressions of the form Equations (1) and (2). Allowing for a variable ∆ρmax
308
results in the more complex dynamics that are captured in our simulations.
309
We also show in Figure S3 the mole fraction of CO2 after 1.5, 2.5, 4.5, 8.5, 11, and
310
13.5 years for binary heterogeneous domains with 16 and 28% volume proportion of high-
311
permeability facies. As a reminder, the ratio between porosities of the two facies is 2.6.
312
This means that grid-cells in the high-permeability facies have 2.6 times more pore volume
313
(PV) than the low-permeability facies grid-cells. Equations (1) and (2) suggest that lower
314
porosity facies have a lower λc and tc . Porosity affects the CO2 transport and the diffusive
315
flux (respectively, Equations (1) and (2) in Text S1), and its dispersion-width (Equation (3)).
316
As shown in Figure S3, the first fingers develop in high-permeability facies (see results of
317
1.5 and 2.5 years suggesting separated island of fingers through high-permeability facies).
318
However, after 4.5 years the lower permeability-porosity facies become saturated with CO2
14
ACS Paragon Plus Environment
Page 14 of 33
Page 15 of 33
Environmental Science & Technology
319
as well. This increases the density contrast and triggers instabilities. As a result, while the
320
advective Darcy flux (∝ k) is lower in low-permeability facies, the finger growth rates (∝ k/φ)
321
can be comparable to that in the higher permeability regions for both cases with 16 and 28%
322
volume proportions of high-permeability facies. This effect becomes more pronounced as the
323
volume proportion of low-permeability facies increases.
324
Figure 2b compares the temporal evolution of σz for different volume proportions of
325
high-permeability facies. The heterogeneous domains with fully connected clusters (28 and
326
40% cases with high connectivity of 27 and 39.6%, respectively) show a very similar evo-
327
lution of σz . The evolution for the 5 and 16% cases exhibit a greater difference because
328
the gravitational flow dynamics are more affected by the spatial organizations of both low-
329
and high-permeability facies. An interesting and novel finding is that the heterogeneous
330
simulations show fairly similar scaling relations as compared to the four flow regimes in
331
their homogeneous counterparts. Specifically, we observe both the first and third diffusion-
332
dominated flow regimes, and that these regimes last the longest for the 5 and 16% cases. The
333
similarities in behavior are clear in Figure 2d, which offers a direct quantitative comparison
334
between the homogeneous and heterogeneous simulations.
335
It is clear that the temporal evolution and the scaling of the third and the fourth regimes
336
are the same for homogeneous and heterogeneous simulations. In the first two flow regimes
337
deviations from the homogeneous case are related to the tortuous pathways arising from
338
the facies architecture. The former has significant practical implications in determining
339
the long-term evolution of CO2 convective mixing. Regardless of the complex fingering
340
and mixing dynamics at early times, the long-term storage permanence is determined by the
341
convective mixing in the fourth flow regime. Based on our findings, if the volume proportions
342
of facies types can be determined and their properties measured (e.g., from well-log data
343
or outcrops), then suitably averaged petrophysical properties (permeability and porosity)
344
should be sufficient to predict the long-term fate and storage performance of CO2 .
345
More generally, consider the variables that determine the mixing dynamics in Equa-
15
ACS Paragon Plus Environment
Environmental Science & Technology
346
tion (1) (and similarly in the Rayleigh number). Because of the low solubility of CO2 in
347
water (∼ 1–3 mol%), the maximum density difference ∆ρmax is only around 1% with modest
348
variation. Brine viscosity depends on temperature, but not significantly on CO2 composi-
349
tion. The diffusion coefficient for CO2 in brine depends on temperature and pressure but
350
also likely varies by less than an order of magnitude. Permeability can vary by many orders
351
of magnitude across different types of formations. However, 1) high permeability aquifers
352
are most favorable for CO2 storage, because they have the highest injectivity, and 2) if the
353
porosity is correlated with permeability at all, the combination φ/k in Equation (1) shows
354
less variation. For these reasons, the on-set time of instability, typical size (wavelength) of
355
fingers, and propagation speed may not vary as much as a dimensionless model may suggest
356
in terms of Rayleigh numbers that vary by many orders of magnitude. Typical Rayleigh
357
numbers for real CO2 storage formations may be more constrained.
358
Below we perform additional sensitivity analyses to better understand the potential effects
359
of heterogeneity and connectivity on the density-driven flow of CO2 . We perform simulations
360
1) with constant average porosities, 2) with within-cell scale anisotropy, 3) with shorter
361
correlation length for facies types, and 4) with permeability variation within each facies.
362
Homogeneous porosity
363
To separate the effects of heterogeneous porosities versus heterogeneous permeabilities, we
364
consider the same heterogeneous permeability distributions as before, but with a homoge-
365
neous (mean, as in Table S1) porosity.
366
Figure 1b shows the mole fraction of CO2 after 1.5, 2.5, 4.5, 8.5, 11, and 13.5 years
367
for a binary domain with 28% high-permeability facies. Media with 1) homogeneous k and
368
homogeneous φ (left column), 2) heterogeneous k and homogeneous φ (middle column), and
369
3) heterogeneous k and heterogeneous φ (right column) are shown. As a linear combination
370
of the porosities of the high- and low-permeability facies (0.34 and 0.13, respectively), the
371
mean porosity values always fall between the two (see Table S1). Therefore a k/φ ratio for
16
ACS Paragon Plus Environment
Page 16 of 33
Page 17 of 33
Environmental Science & Technology
ACS Paragon Plus Environment
Environmental Science & Technology
381
and φ as shown in the left panels. In high-permeability facies, tc and λc are shorter and
382
fingering is more pronounced. In general tc is smaller for all heterogeneous media (regardless
383
of porosity spatial variability) compared to homogeneous ones (Figure 2d).
384
The σz scaling is different for simulations with heterogeneous k but homogeneous φ
385
(Figure 2c). Similar to other cases, σz grows diffusively with σz ∝ t0.5 . The other three
386
distinct flow regimes observed in homogeneous and fully heterogeneous domains are smoothed
387
out in heterogeneous formations with uniform φ. After fingers are triggered the σz grows
388
faster in fully connected clusters, leading to a higher σz for both volume proportions of 28
389
and 40% with fully connected clusters. The volume proportion of 5% shows smaller σz .
390
For volume proportions with fully connected clusters, we see predominantly ballistic flow
391
regime (σz ∝ t1 ). The σz reaches the maximum σz at comparable times as the equivalent
392
homogeneous cases. For 5 and 16% volume proportions (below the percolation threshold), we
393
observe the ballistic regime only in part of the σz scaling. For example, volume proportion of
394
16% has a drop in σz after 10 years. The 5% case also shows some deviations from the ballistic
395
regime. This is related to the spatial organization and the location of high-permeability facies
396
and it is a realization-dependent phenomenon. The ensemble results of many realizations
397
should be considered (the number of realizations are not known a priori and could be as
398
many as hundreds) for cases that are below the percolation threshold (0–20%). However,
399
such analyses are computationally expensive and are outside the scope of this work.
400
Anisotropy ratio
401
As one other potentially important formation property, we consider small-scale (within grid
402
cells) anisotropy. Horizontal permeability, kh , is typically larger than the vertical kv . We
403
define the vertical to horizontal permeability as anisotropy ratio γ = kv /kh , which is typically
404
< 1. We create γ < 1 by increasing the horizontal permeability kh while keeping kv constant.
405
This way we are able to study anisotropy while comparing the results for the same buoyancy
406
scale (i.e., with equal driving force for vertical convective flux) but different γ. There are only
18
ACS Paragon Plus Environment
Page 18 of 33
Page 19 of 33
Environmental Science & Technology
407
a few works that have considered how the dissolution flux of CO2 is affected by anisotropic
408
porous media 68–70 for homogeneous or relatively simple (e.g., layered) heterogeneous media.
409
Figure 3 shows how γ < 1 affects the time evolution of σz in fully heterogeneous binary
410
media, with the 28% volume proportion of high-permeability facies.
411
We find that higher cell-scale anisotropy (i.e., lower γ) leads to earlier onset of convection
412
(i.e., tc ) in both homogeneous and heterogeneous media. However, this effect is more pro-
413
nounced for heterogeneous media. The degree of anisotropy considerably affects the onset
414
time and results in a shorter first diffusive regime. Following the instability onset, we find
415
that γ < 1 does not significantly affect the σz scaling in the early-time convective regime for
416
both homogeneous and heterogeneous media. Anisotropy does appear to lower the σz tem-
417
poral scaling for the last, fully convective, regime (again more pronounced for heterogeneous
418
media). The reason is that, although kv is the same for different γ, kh is higher for smaller
419
γ. This leads to more lateral mixing and diffusion, which lowers the CO2 concentration and
420
density within the sinking fingers. As a result, the (negative) buoyancy driving force and
421
the associated vertical spreading rate are weakened in the last (fourth) regime.
422
Correlation length and variability within facies
423
In the previous sections, we considered an isotropic integral scale (λ) of 10 m (X hereafter).
424
In this section we consider binary structures with 0.5 X and 0.1 X for a volume proportion
425
of 28%. Note that both k and φ are spatially heterogeneous. The results are shown in
426
Figure S4a, and suggest that changes in the correlation length have negligible effect on σz ,
427
and the scaling is similar for all cases.
428
Figure S4b shows σz (for a volume proportion of 28%) with different values of the per-
429
meability variance (σj2 ) within each facies. In previous analyses we did not incorporate any
430
variability in k since we assumed that the spatial organization of facies types is the major
431
source of heterogeneity, and that the spatial variabilities in k and φ within each facies are
432
of secondary importance. 72 Adding σj2 changes our binary systems to bimodal media. Fig-
19
ACS Paragon Plus Environment
Environmental Science & Technology
ACS Paragon Plus Environment
Page 20 of 33
Page 21 of 33
Environmental Science & Technology
441
respectively, 19 and 28% volume proportions of OFC cells are provided in Figure S5. We
442
also provide the results for homogeneous equivalent domains, i.e., with the same mean
443
permeability and porosity, in Figure S6. We note that a mean (scalar) permeability does not
444
capture effective anisotropy that may exist in the complex multiscale heterogeneity (defining
445
such anisotropy requires sophisticated upscaling techniques that are beyond the scope of this
446
work).
447
The difference in fingering behavior between a homogeneous (Figure S6) and a multi-
448
scale heterogeneous domain (Figure S5) appears to be more pronounced than for the binary
449
heterogeneity in terms of onset time and typical wavelength (with more and smaller fingers
450
in the homogeneous domain; see also Figure 4). This would not necessarily be surprising
451
given the wide range and non-trivial distribution of facies permeabilities. However, more
452
interestingly, the global measure of dispersion-width, σz , is surprisingly insensitive to the
453
degree and type of multiscale heterogeneity (Figure 4).
454
The inset in Figure 4 shows the same results but in a linear, rather than logarithmic,
455
scale. It is clear that cases with fully connected clusters (24 and 28%) exhibit identical time
456
evolution for σz . The architecture below the percolation threshold (the 19% case) has a
457
slightly lower growth rate, which reflects how the spatial distribution of disconnected flow
458
paths slows down the sinking of gravitational fingers to the bottom of domain.
459
Clearly, there are also some differences in the σz evolution between the heterogeneous
460
and the ‘equivalent’ homogeneous domains. The flow dynamics in the heterogeneous domain
461
represent essentially a superposition of the flow behavior in different clusters with varying
462
permeabilities. As such, the four flow regimes that are clearly visible in a homogeneous
463
medium are somewhat smeared out in the heterogeneous ones (though still identifiable).
464
A key take-away point is that despite the considerable complexity of the multiscale het-
465
erogeneous media, the scaling behavior for all heterogeneous and homogeneous domains is
466
remarkably similar. Most importantly, once the smallest early-time fingers have percolated
467
through potentially highly tortuous flow paths, temporal evolution of the σz always ends up
21
ACS Paragon Plus Environment
Environmental Science & Technology
ACS Paragon Plus Environment
Page 22 of 33
Page 23 of 33
Environmental Science & Technology
478
years, then it appears that a simple scaling relation is able to provide reasonable estimates,
479
regardless of the complexity of the formation heterogeneity.
480
Acknowledgement
481
The first author was supported by the U. S. Department of Energy (DOE) Office of Fossil
482
Energy funding to Oak Ridge National Laboratory (ORNL) under project FEAA-045. ORNL
483
is managed by UT-Battelle for the U.S. DOE under Contract DE-AC05-00OR22725. The
484
third and fifth authors was supported by funding as part of the Center for Geological Storage
485
of CO2 , an Energy Frontier Research Center funded by DOE, Office of Science, Basic Energy
486
Sciences under Award # DE-SC0C12504.
487
Supporting Information Available
488
The Supporting Information includes the following sections.
489
• Text S1: Governing equations
490
• Table S1: Mean homogeneous permeability and porosity in binary cases
491
• Figure S1: Heterogeneous binary models
492
• Figure S2: Permeability fields of multiscale heterogeneity model
493
• Figure S3: Visualizations of CO2 mole fraction at different simulation times for het-
494
495
496
497
498
erogeneous binary models • Figure S4: Sensitivity results of dispersion-width to correlation length and permeability variance within facies • Figure S5: Visualizations of CO2 mole fraction at different simulation times for multiscale heterogeneous models 23
ACS Paragon Plus Environment
Environmental Science & Technology
499
500
• Figure S6: Visualizations of CO2 mole fraction at different simulation times for homogeneous equivalents of multiscale models
501
• Text S2: References
502
This material is available free of charge via the Internet at http://pubs.acs.org/.
503
References
504
(1) Eccles, J. K.; Pratson, L.; Newell, R. G.; Jackson, R. B. Physical and economic potential
505
of geological CO2 storage in saline aquifers. Environ. Sci. Technol. 2009, 43, 1962–1969.
506
(2) Mathias, S. A.; Gluyas, J. G.; Goldthorpe, W. H.; Mackay, E. J. Impact of maximum
507
allowable cost on CO2 storage capacity in saline formations. Environ. Sci. Technol.
508
2015, 49, 13510–13518.
509
510
(3) Bourg, I. C.; Beckingham, L. E.; DePaolo, D. J. The nanoscale basis of CO2 trapping for geologic storage. Environ. Sci. Technol. 2015, 49, 10265–10284.
511
(4) Yang, C.; Treviño, R. H.; Zhang, T.; Romanak, K. D.; Wallace, K.; Lu, J.; Mickler, P. J.;
512
Hovorka, S. D. Regional assessment of CO2 –solubility trapping potential: A case study
513
of the coastal and offshore Texas Miocene Interval. Environ. Sci. Technol. 2014, 48,
514
8275–8282.
515
516
(5) Orr Jr, F. M. CO2 capture and storage: are we ready? Energy Environ. Sci. 2009, 2, 449–458.
517
(6) Soltanian, M. R.; Amooie, M. A.; Cole, D. R.; Graham, D. E.; Hosseini, S. A.; Hov-
518
orka, S.; Pfiffner, S. M.; Phelps, T. J.; Moortgat, J. Simulating the Cranfield geological
519
carbon sequestration project with high-resolution static models and an accurate equa-
520
tion of state. Int. J. Greenh. Gas Control. 2016, 54, 282–296.
24
ACS Paragon Plus Environment
Page 24 of 33
Page 25 of 33
521
522
Environmental Science & Technology
(7) Juanes, R.; Spiteri, E.; Orr, F.; Blunt, M. Impact of relative permeability hysteresis on geological CO2 storage. Water Resour. Res. 2006, 42, W12418.
523
(8) Gershenzon, N. I.; Ritzi, R. W.; Dominic, D. F.; Soltanian, M.; Mehnert, E.; Ok-
524
wen, R. T. Influence of small-scale fluvial architecture on CO2 trapping processes in
525
deep brine reservoirs. Water Resour. Res. 2015, 51, 8240–8256.
526
(9) Farajzadeh, R.; Salimi, H.; Zitha, P. L.; Bruining, H. Numerical simulation of density-
527
driven natural convection in porous media with application for CO2 injection projects.
528
Int. J. Heat Mass Transfer. 2007, 50, 5054–5064.
529
(10) Rapaka, S.; Chen, S.; Pawar, R. J.; Stauffer, P. H.; Zhang, D. Non-modal growth of
530
perturbations in density-driven convection in porous media. J. Fluid Mech. 2008, 609,
531
285–303.
532
(11) Farajzadeh, R.; Ranganathan, P.; Zitha, P. L. J.; Bruining, J. The effect of heterogeneity
533
on the character of density-driven natural convection of CO2 overlying a brine layer.
534
Adv. Water Resour. 2011, 34, 327–339.
535
536
(12) Hidalgo, J. J.; Dentz, M.; Cabeza, Y.; Carrera, J. Dissolution patterns and mixing dynamics in unstable reactive flow. Geophys. Res. Lett. 2015, 42, 6357–6364.
537
(13) Soltanian, M. R.; Amooie, M. A.; Dai, Z.; Cole, D.; Moortgat, J. Critical dynamics of
538
gravito-convective mixing in geological carbon sequestration. Sci. Rep. 2016, 6, 35921.
539
(14) Agartan, E.; Trevisan, L.; Cihan, A.; Birkholzer, J.; Zhou, Q.; Illangasekare, T. H. Ex-
540
perimental study on effects of geologic heterogeneity in enhancing dissolution trapping
541
of supercritical CO2 . Water Resour. Res. 2015, 51, 1635–1648.
542
543
(15) Daniel, D.; Riaz, A.; Tchelepi, H. A. Onset of natural convection in layered aquifers. J. Fluid Mech. 2015, 767, 763–781.
25
ACS Paragon Plus Environment
Environmental Science & Technology
544
545
(16) Riaz, A.; Hesse, M.; Tchelepi, H.; Orr, F. Onset of convection in a gravitationally unstable diffusive boundary layer in porous media. J. Fluid Mech. 2006, 548, 87–111.
546
(17) Hassanzadeh, H.; Pooladi-Darvish, M.; Keith, D. W. Scaling behavior of convective
547
mixing, with application to geological storage of CO2 . AlChE J. 2007, 53, 1121–1131.
548
(18) Pruess, K. Numerical modeling studies of the dissolution-diffusion-convection process
549
during CO2 storage in saline aquifers. LBNL 2008, 1243E
550
(19) Pau, G. S.; Bell, J. B.; Pruess, K.; Almgren, A. S.; Lijewski, M. J.; Zhang, K. High-
551
resolution simulation and characterization of density-driven flow in CO2 storage in
552
saline aquifers. Adv. Water Resour. 2010, 33, 443–455.
553
(20) Neufeld, J. A.; Hesse, M. A.; Riaz, A.; Hallworth, M. A.; Tchelepi, H. A.; Huppert, H. E.
554
Convective dissolution of carbon dioxide in saline aquifers. Geophys. Res. Lett. 2010,
555
37, L22404.
556
(21) Shahraeeni, E.; Moortgat, J.; Firoozabadi, A. High-resolution finite element methods
557
for 3D simulation of compositionally triggered instabilities in porous media. Comput.
558
Geosci. 2015, 19, 899–920.
559
(22) Soltanian, M. R.; Ritzi, R. W. A new method for analysis of variance of the hydraulic
560
and reactive attributes of aquifers as linked to hierarchical and multiscaled sedimentary
561
architecture. Water Resour. Res. 2014, 50, 9766–9776.
562
(23) Lu, J.; Cook, P.; Hosseini, S.; Yang, C.; Romanak, K.; Zhang, T.; Freifeld, B.;
563
Smyth, R.; Zeng, H.; Hovorka, S. Complex fluid flow revealed by monitoring CO2
564
injection in a fluvial formation. J. Geophys. Res. 2012, 117, B03208.
565
(24) Hosseini, S. A.; Lashgari, H.; Choi, J. W.; Nicot, J.-P.; Lu, J.; Hovorka, S. D. Static and
566
dynamic reservoir modeling for geological CO2 sequestration at Cranfield, Mississippi,
567
USA. Int. J. Greenh. Gas Control. 2013, 18, 449–462. 26
ACS Paragon Plus Environment
Page 26 of 33
Page 27 of 33
Environmental Science & Technology
568
(25) Ritzi, R. W.; Freiburg, J. T.; Webb, N. D. Understanding the (co)variance in petrophys-
569
ical properties of CO2 reservoirs comprising sedimentary architecture. Int. J. Greenh.
570
Gas Control. 2016, 51, 423–434.
571
(26) Gershenzon, N. I.; Soltanian, M. R.; Ritzi, R. W.; Dominic, D. F. Influence of small
572
scale heterogeneity on CO2 trapping processes in deep saline aquifers. Energy Procedia
573
2014, 59, 166–173.
574
(27) Krevor, S.; Pini, R.; Zuo, L.; Benson, S. M. Relative permeability and trapping of CO2
575
and water in sandstone rocks at reservoir conditions. Water Resour. Res, 2012, 48,
576
W02532.
577
(28) Dance, T.; Paterson, L. Observations of carbon dioxide saturation distribution and
578
residual trapping using core analysis and repeat pulsed-neutron logging at the CO2
579
CRC Otway site. Int. J. Greenh. Gas Control. 2016, 47, 210–220.
580
(29) Dai, Z.; Viswanathan, H.; Middleton, R.; Pan, F.; Ampomah,W.; Yang, C.; Jia, W.;
581
Xiao, T.; Lee, S.-Y.; McPherson, B. CO2 accounting and risk analysis for CO2 seques-
582
tration at enhanced oil recovery sites. Environ. Sci. Technol. 2016, 50, 7546–7554.
583
(30) Trevisan, L.; Pini, R.; Cihan, A.; Birkholzer, J. T.; Zhou, Q.; Illangasekare, T. H.
584
Experimental analysis of spatial correlation effects on capillary trapping of supercritical
585
CO2 at the intermediate laboratory scale in heterogeneous porous media. Water Resour.
586
Res. 2015, 51, 8791–8805.
587
(31) Trevisan, L.; Pini, R.; Cihan, A.; Birkholzer, J. T.; Zhou, Q.; González-Nicolás, A.;
588
Illangasekare, T. H. Imaging and quantification of spreading and trapping of carbon
589
dioxide in saline aquifers using meter-scale laboratory experiments. Water Resour. Res.
590
2016, 53, 485–502
591
(32) Trevisan, L.; Krishnamurthy, P.; Meckel, T. Impact of 3D capillary heterogeneity and
27
ACS Paragon Plus Environment
Environmental Science & Technology
592
bedform architecture at the sub-meter scale on CO2 saturation for buoyant flow in
593
clastic aquifers. Int. J. Greenh. Gas Control. 2017, 56, 237–249.
594
595
(33) Stauffer, D.; Aharony, A. Introduction to Percolation Theory; Taylor and Francis, Bristol, Pa., 1992.
596
(34) Sahimi, M. Applications of Percolation Theory; Taylor and Francis, Bristol, Pa., 1994.
597
(35) Guin, A.; Ritzi, R. W. Studying the effect of correlation and finite-domain size on
598
spatial continuity of permeable sediments. Geophys. Res. Lett. 2008, 35, L10402.
599
(36) Proce, C. J.; Ritzi, R. W.; Dominic, D. F.; Dai, Z. Modeling multiscale heterogeneity
600
and aquifer interconnectivity. Groundwater 2004, 42, 658–670.
601
(37) Ramanathan, R.; Guin, A.; Ritzi, R. W.; Dominic, D. F.; Freedman, V. L.;
602
Scheibe, T. D.; Lunt, I. A. Simulating the heterogeneity in braided channel belt de-
603
posits: 1. A geometric-based methodology and code. Water Resour. Res. 2010, 46,
604
W04516.
605
(38) Guin, A.; Ramanathan, R.; Ritzi, R. W.; Dominic, D. F.; Lunt, I. A.; Scheibe, T. D.;
606
Freedman, V. L. Simulating the heterogeneity in braided channel belt deposits: 2.
607
Examples of results and comparison to natural deposits. Water Resour. Res. 2010, 46,
608
W04515.
609
610
611
612
(39) Knudby, C.; Carrera, J. On the relationship between indicators of geostatistical, flow and transport connectivity. Adv. Water Resour. 2005, 28, 405–421. (40) Silliman, S. The influence of grid discretization on the percolation probability within discrete random fields. J. of Hydrol. 1990, 113, 177–191.
613
(41) Gershenzon, N.; Soltanian, M.; Ritzi, R. W.; Dominic, D. F. Understanding the impact
614
of open-framework conglomerates on water–oil displacements: the Victor interval of the
615
Ivishak Reservoir, Prudhoe Bay Field, Alaska. Pet. Geosci. 2015, 21, 43–54. 28
ACS Paragon Plus Environment
Page 28 of 33
Page 29 of 33
Environmental Science & Technology
616
(42) Massabó, M.; Bellin, A.; Valocchi, A. Spatial moments analysis of kinetically sorbing
617
solutes in aquifer with bimodal permeability distribution. Water Resour. Res. 2008,
618
44, W09424.
619
620
621
622
(43) Carle, S. F. T-PROGS: Transition probability geostatistical software. Univ. of Calif., Davis, Calif. 1999. (44) Bernabé, Y.; Mok, U.; Evans, B. in Thermo-Hydro-Mechanical Coupling in Fractured Rock, Birkhäuser Basel, 2003, 937–960.
623
(45) Deng, H.; Stauffer, P. H.; Dai, Z.; Jiao, Z.; Surdam, R. C. Simulation of industrial-scale
624
CO2 storage: Multi-scale heterogeneity and its impacts on storage capacity, injectivity
625
and leakage. Int. J. Greenh. Gas Control. 2012, 10, 397–418.
626
627
628
629
(46) Huang, L.; Ritzi, R. W.; Ramanathan, R. Conservative models: Parametric entropy vs. temporal entropy in outcomes. Groundwater 2012, 50, 199–206. (47) Lunt, I.; Bridge, J.; Tye, R. A quantitative, three-dimensional depositional model of gravelly braided rivers. Sedimentology 2004, 51, 377–414.
630
(48) Gershenzon, N. I.; Soltanian, M. R.; Ritzi, R. W.; Dominic, D. F.; Keefer, D.; Shaf-
631
fer, E.; Storsved, B. How does the connectivity of open-framework conglomerates within
632
multi-scale hierarchical fluvial architecture affect oil-sweep efficiency in waterflooding?
633
Geosp. 2015, 11, 2049–2066.
634
(49) Moortgat, J.; Sun, S.; Firoozabadi, A. Compositional modeling of three-phase flow
635
with gravity using higher-order finite element methods. Water Resour. Res. 2011, 47,
636
W05511.
637
638
(50) Moortgat, J.; Firoozabadi, A. Higher-order compositional modeling with Fickian diffusion in unstructured and anisotropic media. Adv. Water Resour. 2010, 33, 951–968.
29
ACS Paragon Plus Environment
Environmental Science & Technology
639
(51) Moortgat, J.; Li, Z.; Firoozabadi, A. Three-phase compositional modeling of CO2 in-
640
jection by higher-order finite element methods with CPA equation of state for aqueous
641
phase. Water Resour. Res. 2012, 48, W12511.
642
(52) Moortgat, J.; Firoozabadi, A. Mixed-hybrid and vertex-discontinuous-Galerkin finite el-
643
ement modeling of multiphase compositional flow on 3D unstructured grids. J. Comput.
644
Phys. 2016, 315, 476–500.
645
(53) Moortgat, J. B.; Firoozabadi, A.; Li, Z.; Espósito, R. CO2 injection in vertical and
646
horizontal cores: Measurements and numerical simulation. SPE J. 2013, 18, 331–344.
647
(54) Moortgat, J.; Amooie, M. A.; Soltanian, M. R. Implicit finite volume and discontinuous
648
Galerkin methods for multicomponent flow in unstructured 3D fractured porous media.
649
Adv. Water Resour. 2016, 96, 389–404.
650
651
(55) Amooie, M. A.; Soltanian, M. R.; Moortgat, J. Hydro-thermodynamic mixing of fluids across phases in porous media. Geophys. Res. Lett. 2017, 44,3624–3634.
652
(56) Amooie, M. A.; Soltanian, M. R.;Xiong, F.; Dai, Z.; Moortgat, J. Mixing and spread-
653
ing of multiphase fluids in heterogeneous bimodal porous media. Geomechanics and
654
Geophysics for Geo-Energy and Geo-Resources 2017, 1–20.
655
656
657
658
(57) Moortgat, J. Viscous and gravitational fingering in multiphase compositional and compressible flow. Adv. Water Resour. 2016, 89, 53–66 (58) Firoozabadi, A. Thermodynamics and Applications of Hydrocarbon Energy Production; McGraw-Hill Education, New York, 2016.
659
(59) Nasrabadi, H.; Moortgat, J.; Firoozabadi, A. New three-phase multicomponent com-
660
positional model for asphaltene precipitation during CO2 injection using CPA-EOS.
661
Energy Fuels 2016, 30, 3306–3319.
30
ACS Paragon Plus Environment
Page 30 of 33
Page 31 of 33
662
663
664
665
666
667
Environmental Science & Technology
(60) Li, Z.; Firoozabadi, A. Cubic-plus-association equation of state for water-containing mixtures: Is cross association necessary? AIChE J. 2009, 55, 1803–1813. (61) Li, Z.; Firoozabadi, A. General strategy for stability testing and phase-split calculation in two and three phases. SPE J. 2012, 17, 1–096. (62) Hidalgo, J. J.; Carrera, J.; Medina, A. Role of salt sources in density-dependent flow. Water Resour. Res. 2009, 45, W05503.
668
(63) Emami-Meybodi, H.; Hassanzadeh, H.; Green, C. P.; Ennis-King, J. Convective disso-
669
lution of CO2 in saline aquifers: Progress in modeling and experiments. Int. J. Greenh.
670
Gas Control. 2015, 40, 238–266.
671
672
673
674
675
676
(64) Martinez, M. J.; Hesse, M. A. Two-phase convective CO2 dissolution in saline aquifers. Water Resour. Res. 2016, 52, 585–599. (65) Aris, R. On the dispersion of a solute in a fluid flowing through a tube. Proc. R. Soc. Lond. A 1956, 235, 67–77. (66) Le Borgne, T.; Dentz, M.; Carrera, J. Lagrangian statistical model for transport in highly heterogeneous velocity fields. Phys. Rev. Lett. 2008, 101, 090601.
677
(67) Buscheck, T. A.; White, J. A.; Carroll, S. A.; Bielicki, J. M.; Aines, R. D. Managing
678
geologic CO2 storage with pre-injection brine production: a strategy evaluated with a
679
model of CO2 injection at Snøhvit. Energy Environ. Sci. 2016, 9, 1504–1512.
680
(68) Cheng, P.; Bestehorn, M.; Firoozabadi, A.Effect of permeability anisotropy on
681
buoyancy-driven flow for CO2 sequestration in saline aquifers.Water Resour. Res. 2012,
682
48, .
683
684
(69) Green, C.; Ennis-King, J. Steady dissolution rate due to convective mixing in anisotropic porous mediaAdv. Water Resour. 2014, 73, 65–73.
31
ACS Paragon Plus Environment
Environmental Science & Technology
685
686
687
688
689
690
(70) De Paoli, M.; Zonta, F.; Soldati, A. Dissolution in anisotropic porous media: Modelling convection regimes from onset to shutdown. Physics of Fluids 2017, 29, . (71) Divins, D. Total Sediment Thickness of the World’s Oceans & Marginal Seas; NOAA National Geophysical Data Center, Boulder, CO, 2003 (72) Desbarats, A. J. Numerical estimation of effective permeability in sand-shale formations. Water Resour. Res. 1987, 23, 273–286.
32
ACS Paragon Plus Environment
Page 32 of 33
Page 33 of 33
Environmental Science & Technology
ACS Paragon Plus Environment