Currently, the technology of CO
_{2}
capture and storage (CCS) has become the main issue for climate change and global warming. Among CCS technologies, the prediction of CO
_{2}
behavior underground is very critical for CO
_{2}
storage design, especially for its safety. Hence, the purpose of this paper is to model and simulate CO
_{2}
flow and its heat transfer characteristics in a storage site, for more accurate evaluation of the safety for CO
_{2}
storage process. In the present study, as part of the storage design, a micro porescale model was developed to mimic real porous structure, and computational fluid dynamics was applied to calculate the CO
_{2}
flow and thermal fields in the micro porescale porous structure. Three different configurations of 3dimensional (3D) micropore structures were developed, and compared. In particular, the technique of assigning random pore size in 3D porous media was considered. For the computation, physical conditions such as temperature and pressure were set up, equivalent to the underground condition at which the CO
_{2}
fluid was injected. From the results, the characteristics of the flow and thermal fields of CO
_{2}
were scrutinized, and the influence of the configuration of the micropore structure on the flow and scalar transport was investigated.
1. Introduction
To overcome global warming, considerable international attention is being paid these days to the technology of CO
_{2}
capture and storage (CCS)
[1

4]
. But CO
_{2}
leakage, or the fracture of the geological formation where CO
_{2}
is stored, presents a serious problem. Hence, for the commercialization of CCS, one of the most important things is to obtain an accurate prediction or simulation of the behavior of CO
_{2}
in its geological storage layer. For these purposes, numerical simulation methods can be one of the solutions for predicting the multiphase flow and thermal characteristics of CO
_{2}
, near the critical temperature and pressure conditions of underground. To date, simple numerical methods using Darcy’s law have been widely used for simulation techniques
[5]
. On the other hand, the lattice Boltzmann (LB) method has been developed in the field of earth science for molecular level simulation
[6
,
7]
. However, the application of computational fluid dynamics (CFD) is very rare, in spite of its moderate applicability for predicting CO
_{2}
behavior with increased accuracy and low computational cost, compared with Darcy’s law and LB methods, respectively. Hence, in the present study, CFD is applied for calculating CO
_{2}
behavior in a microporous media, which represents an underground storage layer. The present CFD results can be applied to the design of CCS process, especially for CO
_{2}
injection processes, providing permeability prediction with high fidelity for a storage site, etc., which is the final goal of the present study.
For this purpose, a micro porescale model is developed, and this consists of poreholes, which are surrounded by grains. Using this model, a 3dimensional (3D) computational grid is generated and tested for the simulation of CO
_{2}
flow, changing the micro porehole size and its distribution, to investigate the effect of pore configuration on the flow and thermal fields. In particular, to replicate the real porous structure of a storage layer, the porehole size and distribution are randomly given. Near critical pressure and temperature conditions for CO
_{2}
, the evolution of the CO
_{2}
flow and thermal fields are investigated by varying the micropore structure. This will be very helpful for designing optimal CO
_{2}
injection and sequestration systems.
2. Materials and Methods
 2.1. Numerical Procedure
 2.1.1. Governing equations
The continuity and momentum equations for incompressible fluid are expressed as follows:
The computational domain and grid allocation for the three different cases.
Here,
u_{i}
,
p
and τ
_{ij}
are velocity, pressure and viscous stress tensor, respectively. The viscous stress tensor τ
_{ij}
is defined as τ
_{ij}
= 2
μS_{ij}
– (2/3)
μ
(∂
u_{k}
/ ∂
х_{k}
)
δ_{ij}
, and
S_{ij}
is the rate of strain tensor, defined as
S_{ij}
= 0.5 (∂
u_{i}
/ ∂
х_{j}
+ ∂
u_{j}
/ ∂
х_{i}
).
The energy equation to represent the evolution of a passive scalar is expressed as follows:
Where,
T
is temperature,
c
is fluid specific heat and
k
is fluid conductivity. For the spatial discretization of Eqs. (13), the secondorder UPWIND scheme is used. To avoid pressurevelocity decoupling, a SIMPLE algorithm is applied. In the present study, an unstructured grid system is adopted using polyhedral meshes, and STARCCM+ ver. 3.04 is used to solve Eqs. (13)
[8]
. Before conducting the main calculation, flow through a lattice flow cell model
[9]
was computed, to evaluate the present numerical methods. Comparing our results with that of Mazaheri et al.
[9]
, computational results indirectly show the validity of the present calculation
[10]
.
 2.1.2. Numerical methods
Fig. 1
shows a part of the computational domain for various configurations to investigate the effect of the porous media configuration on flow and thermal fields. For reference, a part of the
Characteristics of the three different computational domains
Characteristics of the three different computational domains
computational domain is represented in the figure, to show pore configuration clearly. As shown in Choi et al.
[10]
, a pore structure exists among grains. Hence, a pore represents the fluid flow passage surrounded by grains, and CO
_{2}
flow passes through the micropores in the present study. To compare the effect of micropore size, shape and distribution, three different cases of pore structure are selected, as described in
Table 1
.
For reference, the grain size in
Table 1
is determined from the wellknown experiment
[11]
, which found out the permeability of sandstones. For porosity, the value of 0.20.3 is usually used for
The construction of the computational domain for case 3. (a) Random distribution of various size grains, (b) pore structure after subtracting grains from the cubic volume.
Permeability of the three different cases compared with the experimental results [16].
the modeling of aquifer
[12

15]
where CO
_{2}
may be sequestrated, hence the porosity is adopted around the range in the present study, considering experiment
[11]
as well. Also, the permeability of the present 3D model is compared with experiment
[11]
in
Fig. 3
, and this will be discussed in section 3 in detail.
Cases 1 and 2 have uniform grain sizes and distribution, respectively. It is noted that random packing and the irregular size of the grains approaches the more realistic situation of an underground porous structure similar to rock or sandstone. Hence, case 3 is further developed for a complicated model considering anisotropy of the underground pore structure. For this purpose, a modeling technique is developed and tested
[16]
, which gives a random distribution of grain size and location. Case 3 has a random distribution of four different grain sizes, to mimic the real underground situation.
Fig. 2
shows how to construct the size distribution of grains for case 3, resulting in anisotropic micropore structure. As depicted in
Fig. 2
(a) with different colors, the random distribution of grains having four different grain diameters is applied in a hexahedral cube. After that, the grains are subtracted from the volume of the cube, and then the micropore structure remains. This micropore structure is used for the computational domain for case 3.
Contours of streamwise velocity for the three different cases.
For the flow boundary conditions, the CO
_{2}
flow enters through the inlet, which is located at the left plane perpendicular to the xaxis of the cubic array of micropores, and comes out from the outlet (perpendicular to the right plane), to simulate CO
_{2}
movement through the underground. For walls, the noslip condition is used, and a symmetry condition is applied to the spanwise and transverse directions. For the thermal boundary condition, an isothermal condition is adopted for the walls, which are heated by higher temperature than that of the fluid. For the grid allocation, unstructured grids with polyhedral meshes are applied, and the grid number is carefully selected around 1,000,000. For reference, to consider the supercritical condition of CO
_{2}
in an underground position, pressure is set to 100 bar and the temperature conditions of inlet and wall are given as 313 and 332 K, respectively. The variation of physical properties of CO
_{2}
is considered, such as density, specific heat, conductivity, and viscosity of CO
_{2}
. It is noted that the change of the physical properties is very noticeable, even if the temperature is changed in the order of 10 K. This may affect the CO
_{2}
flow and final thermal fields, and has to be considered.
3. Results and Discussion
Before conducting the main calculation, the permeability of sandstone is computed to validate our modeling of micropore structure, using the three different cases in
Table 1
. The permeability is calculated using the following equation as
κ
=
νμ
(
dх
/
dP
). Here,
κ
is permeability,
ν
is superficial fluid flow velocity and
μ
is kinematic viscosity of fluid.
Fig. 3
shows the permeability distribution of sandstones that have different porosity, and our data compared with the experimental results
[11]
. As can be seen in the figure, the permeability of case 3 is very close to the experimental one; however the results of cases 1 and 2 are higher than that of the experimental case. Hence, it is noted that the random distribution of four different grain sizes may be more realistic for the underground situation. We now proceed to discuss the effect of the configuration of micropore structure on the flow and thermal fields.
Figs. 4
and
5
show the contours of streamwise velocity and skin friction coefficient for a crosssectional plane slicing the
Contours of skin friction coefficient for the three different cases.
Distribution of crosssectional area along the streamwise direction.
Contours of temperature for the three different cases.
middle of the cube, respectively. In case 3, the magnitude of the streamwise velocity is higher than those of the other cases. Furthermore, the location of higher velocity magnitude is distributed irregularly in case 3, but cases 1 and 2 have uniform higher locations, where the gap between grains becomes narrower. It is interesting that case 1 consists of the smallest grains among three cases over the entire computational domain, but the highest velocity magnitude appears in case 3. This can be explained by the crosssectional area along the streamwise direction, as in
Fig. 6
. Case 3 has the higher amplitude of area variation, compared with the others.
Accordingly, the magnitude of the skin friction coefficient is higher where the flow is accelerated in
Fig. 5
. The skin friction coefficients of cases 1 and 2 have a higher value at the throat region, where the velocity is increased between grains, and are distributed regularly. In case 3, the higher magnitude of skin friction coefficient is also located at the throat region; however the throat region is randomly distributed, resulting in irregular dis
Contours of wall Nusselt number for the three different cases.
tribution of the higher magnitude coefficients. This may affect the permeability, as illustrated in
Fig. 3
, and finally have a great influence on the evolution of the scalar field, as will be discussed in the following.
Figs. 7
and
8
show the contours of temperature and wall Nusselt number for the three different cases. For all the cases, the temperature is higher where the magnitude of the streamwise velocity is lower. In particular, approaching the outlet, the temperature is increased in the order from case 2, via case 1, to case 3. To look into these phenomena, the contour of wall Nusselt number is illustrated in
Fig. 8
. Comparing the three different cases, the magnitude of the wall Nusselt number is larger, and its distribution is irregular in case 3, because of the complex flow pattern. This can be confirmed from the following figure.
Fig. 9
shows the isosurface of Λ
_{2}
value and contour of wall Nusselt number for the three different cases. To define the vortices in the present study, a negative Λ
_{2}
method is used to capture a vortical flow region, as proposed by Jeong and Hussain
[17]
. So, at every grid point calculation was made for finding a quantity, Λ
_{2}
, the second largest one among the three eigenvalues of ? (1/
ρ
)(∂
^{2}
P
/ ∂
х_{i}
∂
х_{i}
), or specially of its equivalent
S_{ik}
S_{kj}
+ Ω
_{ik}
Ω
_{kj}
, where Ω
_{ik}
and Ω
_{kj}
are the vorticity tensor, and
S_{ik}
and
S_{kj}
the strain rate tensor. attention is paid to the value of Λ
_{2}
calculated at grid points over the whole computational domain, and the negative Λ
_{2}
region is regarded as a vortical flow region. For cases 1 and 2, the vortical structure and the higher region of wall Nusselt number appear regularly, and this regular flow pattern induces a uniform heat transfer pattern. Also, the higher region is confined to a small region. However, the generation of the vortical structure is irregular, and the corresponding wall Nusselt number is randomly distributed in case 3. Furthermore, the higher region of wall Nusselt Number is more broadly distributed, compared with the other cases. This can mimic the flow and thermal transport phenomena in a real porestructure. Hence, in the present study, it can be found that the modeling of porestructure can greatly affect the prediction of the flow and thermal fields. This is very important in developing an accurate micropore model for the application of the optimal design of CO
_{2}
injection and storage.
Isosurface of Λ_{2} and contour of wall Nusselt number for the three different cases.
4. Conclusions
In the present study, CFD is applied for calculating CO
_{2}
behavior in a microporous media, which represents a storage layer inside underground. In particular, a micro porescale model is developed, and this consists of poreholes that are surrounded by grains. To investigate the effect of the porestructure configuration, three different allocations of micrograins are considered. In particular, the porehole size and distribution are randomly given, to replicate the real porous structure of the storage layer. For the same size and regular arrangement of grains, the characteristics of the surface friction and heat transfer show a similar repeating pattern. However, if the micrograins are randomly allocated, complex flow and thermal fields appear. It can be found that the modeling of pore structure can greatly affect the prediction of its flow and thermal fields. This is very important in developing a micropore model for the application of the optimal design of CO
_{2}
injection and storage in the CCS process. For example, using the CFD technique, the permeability of a micropore structure can be calculated, which mimics the real porous structure in the storage site. For reference, this numerical prediction is much better than that of the simple calculation by Darcy’s law
[16]
, which is still generally used in the design and engineering of the CCS process. Then, the predicted permeability with high fidelity can be applied to the system design of CCS process, in particular for CO
_{2}
injection design, and this result in a decrease of cost, and increase of safety.
 Nomenclature

c fluid specific heat

kfluid thermal conductivity

Nu wall Nusselt number

p pressure

Sijthe rate of strain tensor

T fluid temperature

uivelocity

κ permeability

μ kinematic molecular viscosity

τijviscous stress tensor

ρ fluid density
Acknowledgements
This work was conducted on behalf of the Ministry of Land, Transport and Maritime Affairs (MLTM) of Korean government under their “Development of Technology for CO2Ocean Sequestration” program. The authors would like to thank all members of our collaborative research project.
Gaspar Ravagnani AT
,
Ligero EL
,
Suslick SB.
2009
CO2sequestration through enhanced oil recovery in a mature oil field.
J. Pet. Sci. Eng.
65
129 
138
Bachu S.
2002
Seqeustration of CO2in Geological media in response to climate change: road map for site selection using the transform of the geological space into the CO2phase space.
Energy Convers. Manag.
43
87 
102
Borchiellini R
,
Massardo AF
,
Santarelli M.
2002
Carbon tax vs CO2sequestration effects on environomic analysis of existing power plants.
Energy Convers. Manag.
43
1425 
1443
Keith D
,
Lavoie R.
2009
An overview of the wabamun area CO2sequestration project (WASP).
Energy Procedia
1
2817 
2824
Helmig R
,
Bastian P
,
Class H
1998
Architecture of the modular program system MUFTEUG for simulating multiphase flow and transport processes in heterogeneous porous media.
Mathemat. Geol.
2
123 
131
Kang Q
,
Tsimpanogiannis IN
,
Ahang D
,
Lichtner PC.
2005
Numerical modeling of porescale phenomena during CO2sequestration in oceanic sediments.
Fuel Process. Technol.
86
1647 
1665
Suekane T
,
Soukawa S
,
Iwatani S
,
Tsushima S
,
Hirai S.
2005
Behavior of supercritical CO2injected into porous media containing water.
Energy
30
2370 
2382
2006
STARCCM+ ver. 3.0 user guide.
CDadapco
Melville
Mazaheri AR
,
Zerai B
,
Ahmadi G
2005
Computer simulation of flow through a lattice flow  cell model.
Adv. Water Resour.
28
1267 
1279
Choi HS
,
Choi YS
,
Park HC
2010
The characteristics of CO2flow and thermal field in a porous media.
American Society of Mechanical Engineers
Washington
Proceedings of the 14th International Heat Transfer Conference (IHTC14)
Washington, DC.
2010 Aug 813
Paper no. IHTC1423365
983 
988
Adler PM
,
Jacquin CG
,
Quiblier JA.
1990
Flow in simulated porous media.
Int. J. Multiph. Flow
16
691 
712
Mosthaf K.
2007
CO2 strage into dipped saline aquifers including ambient water flow [dissertation].
Universitat Stuttgart
Stuttgart
Li Q
,
Wu Z
,
Li X.
2009
Prediction of CO2leakage during sequestration into marine sedimentary strata.
Energy Convers. Manag.
50
503 
509
Birkholzer JT
,
Zhou Q
,
Tsang CF.
2009
Largescale impact of CO2storage in deep saline aquifers: a sensitivity study on pressure response in stratified systems.
Int. J. Greenh. Gas Control
3
181 
194
Wu YS
,
Pan L
,
Pruess K.
2004
A physically based approach for modeling multiphase fracture?matrix interaction in fractured porous media.
Adv. Water Resour.
27
875 
887
Choi HS.
2011
Numerical simulation of supercritical CO2flow in a geological storage reservoir of ocean.
J. Korean Soc. Environ. Eng.
34
251 
257
Jeong J
,
Hussain F.
1995
On the identification of a vortex.
J. Fluid Mech.
285
69 
94