Article pubs.acs.org/EF
Numerical Predictions of Experimentally Observed Methane Hydrate Dissociation and Reformation in Sandstone Knut A. Birkedal,*,∥,⊥,† C. Matt Freeman,⊥,‡ George J. Moridis,⊥,§ and Arne Graue⊥,† †
Department of Physics and Technology, University of Bergen, 5200 Bergen, Norway Hilcorp Energy Company, Houston, Texas 77002, United States § Lawrence Berkeley National Laboratory, Berkeley, California 94720, United States ‡
S Supporting Information *
ABSTRACT: Numerical tools are essential for the prediction and evaluation of conventional hydrocarbon reservoir performance. Gas hydrates represent a vast natural resource with a significant energy potential. The numerical codes/tools describing processes involved during the dissociation (induced by several methods) for gas production from hydrates are powerful, but they need validation by comparison to empirical data to instill confidence in their predictions. In this study, we successfully reproduce experimental data of hydrate dissociation using the TOUGH+HYDRATE (T+H) code. Methane (CH4) hydrate growth and dissociation in partially water- and gas-saturated Bentheim sandstone were spatially resolved using Magnetic Resonance Imaging (MRI), which allows the in situ monitoring of saturation and phase transitions. All the CH4 that had been initially converted to gas hydrate was recovered during depressurization. The physical system was reproduced numerically, using both a simplified 2D model and a 3D grid involving complex Voronoi elements. We modeled dissociation using both the equilibrium and the kinetic reaction options in T+H, and we used a range of kinetic parameters for sensitivity analysis and curve fitting. We successfully reproduced the experimental results, which confirmed the empirical data that demonstrated that heat transport was the limiting factor during dissociation. Dissociation was more sensitive to kinetic parameters than anticipated, which indicates that kinetic limitations may be important in short-term core studies and a necessity in such simulations. This is the first time T+H has been used to predict empirical nonmonotonic dissociation behavior, where hydrate dissociation and reformation occurred as parallel events.
■
tests.20,11,21,22 A series of studies have also aimed at codevalidation through comparison with experimental data.23−31 This study compares experimental dissociation data with numerical data in order to determine the robustness and advantages of numerical modeling. T+H32 is a multicomponent, multiphase fluid and heat flow simulator that describes hydrate formation and dissociation through coupled equations of mass and heat balance. These processes are predicted by an equilibrium and kinetic model. Phase transitions in the equilibrium model are governed by temperature and pressure, where the different phase combinations are described in Figure 1.7 The kinetic model is described by an Arrhenius-type equation, where kinetic parameters are considered in addition to temperature and pressure. In a study involving several gas hydrate configurations and long-term production, the equilibrium and the kinetic models exhibited a similar dissociation response.33 In short-term processes (such as laboratory studies of dissociation in core samples that lasted a few hours), the kinetic and the equilibrium models exhibited significantly different behavior, which converged as time advanced. Thus, Kowalsky and Moridis (2007)33 recommended the use of the kinetic model in short-term processes. The kinetic model is based on the work of Bishnoi and coworkers34,35 but was modified to account for geometry
INTRODUCTION Gas hydrate is a solid state of special gases (capable of forming hydrates) and water, in which hydrogen bonding and van der Waals interaction forces provide stability to the crystal structure. Of the various hydrate-forming gases in natural hydrates, methane (CH4) is overwhelmingly the most common. The CH4-hydrate stability is governed by pressure and temperature, and its hydration reaction is commonly described by CH4 + NH H 2O = CH4 ·NH H 2O + heat
where the hydration number NH is estimated to be 5.99 ± 0.07.1 Vast resources are associated with gas hydrates,2 and depressurization is often considered one of the more promising approaches for hydrate dissociation.3 Dissociation requires less energy relative to thermal methods but is limited by the endothermic reaction and therefore relies on heat transfer from the surroundings. Numerical modeling is a valuable tool for conventional oil and gas production. Space and time are discretized into a finite number of subdivisions, where coupled equations describe the physical and chemical conditions of each subdomain. A range of simulation tools are available for hydrate simulation.4 Several case studies have been performed due to limited field data,5−10 while others use field specific parameters for evaluation.11−19 Code verification has been performed through comparison with available production data and open-hole pressure-response © 2014 American Chemical Society
Received: January 21, 2014 Revised: August 20, 2014 Published: August 25, 2014 5573
dx.doi.org/10.1021/ef500255y | Energy Fuels 2014, 28, 5573−5586
Energy & Fuels
Article
consequently, the reaction surface area. Less cooling is required in the presence of nucleation surfaces, as the contact-angle reduces the Gibbs free energy and increases the probability of growth. Minerals may therefore trigger growth as preferred nucleation sites. Capillary pressure and fluid distribution within the pore space is continuously changing during hydrate growth and dissociation, and the hydrate crystal may also act as nucleation site itself. Thus, the kinetic parameters are expected to have substantially different values in a nonagitated porous environment. Some of these mechanisms are saturationdependent, and a time-variable kinetic approach may be required to fully capture the complexity of gas hydrate formation and dissociation. Fundamental understanding of processes in hydrate reservoirs during dissociation is imperative prior to full-scale longterm gas hydrate production. Building confidence based on good agreement between experimental data and numerical predictions is therefore essential in order to bridge knowledge gaps and improve our understanding of the complex behavior of gas hydrates. The main purpose of this study is to address this issue by comparing the numerical predictions of the T+H code with laboratory-measured and obtained dissociation data, augmented by visual observations from in situ MR imaging that has been used to monitor the processes occurring within the porous media during the experiments.
Figure 1. Hydrate pressure and temperature equilibrium diagram for different phase combinations.
irregularities in the assumed grain sphericity. The final mass accumulation term is given by QH =
■
∂M = −K 0eΔEa / RT FAA(feq − fv ) ∂t
EXPERIMENTAL DESCRIPTION
Experimental Setup. The experimental setup (Figure 2) involved a Varian Unity/Inova 85.7 MHz MRI that could spatially resolve water saturation and subsequently hydrate growth patterns in a porous medium sample. Hydrogen in water and CH4 were detected through a standard spin−echo sequence. The resonance of hydrogen in the highly organized hydrate crystal is short-lived and is not detectable through the spin echo sequence. Thus, changes in hydrate saturation were given implicitly through loss of signal. We used two Quizix C6000-10K-HC-HT pumps for line (CH4) pressure and permeability measurements, where three Paroscientific Digiquartz pressure transducers logged absolute and differential pressure. The confining fluid was Fluorinert FC-40, selected because it contains no hydrogen atoms and it minimizes RF loss because of its low dielectric properties. A Quizix QX-6000 pump pressurized the Fluorinert FC-40 fluid, and a reciprocating pump circulated the confining fluid in the high-pressure system. The confining tubing also involved PVC lines in which
where K0 is the intrinsic hydration reaction constant [kg/m2·Pa· s], ΔEa is the hydration activation energy [J/mol], R is the universal gas constant [8.314 J/(molK)], T is temperature [K], FA is an area adjustment factor to correct for irregularities in grain shape and differences in hydrate shape and distribution,32 A is the active reaction surface area [m2], feq is the fugacity [Pa] at equilibrium at temperature T, and f v is the vapor fugacity [Pa] at temperature T. The activation energy and intrinsic rate was measured by Bishnoi and co-workers34,35 for a pure hydrate system using stirred tank semibatch reactors. Growth in the presence of solid surfaces is assumed to be substantially different, both in terms of kinetics and stability.36 The pore size distribution and the mineral geometry will affect the fluid distribution, and,
Figure 2. Experimental setup used during this study (Figure modified with permission from ref 37. Copyright 2008, University of Bergen.). 5574
dx.doi.org/10.1021/ef500255y | Energy Fuels 2014, 28, 5573−5586
Energy & Fuels
Article
Figure 3. Cross-sectional view of the composite core holder. Image modified with permission from ref 37. Copyright 2008, University of Bergen. antifreeze circulated for temperature control of the Fluorinert. The confining temperature was measured inside the core holder (Figure 3) using a type-T thermocouple with 0.1 °C resolution, and was maintained at 4 ± 0.3 °C during dissociation. Temperature measurements within the sample were not possible for this experimental setup. The custom-built composite fiberglass high-pressure cell we used for this study was compatible with the magnetic field. The core holder and other properties of the experimental setup have been described elsewhere37 along with some results from one of the depressurization tests.38 They will not be discussed in detail here, but are synoptically illustrated in Figure 3. We also used a second experimental setup without imaging capabilities for test 3. In this system, the pressure vessel was a RHC-Hassler type core holder with aluminum housing, Buna-N sleeve, stainless steel end-pieces and distribution plugs. This setup also involved high precision pumps and accurate pressure transducers (see Figure I, Supporting Information for further description of this experimental setup). Water mixed with antifreeze was used as cooling fluid and temperature was measured both in the confining fluid and at the core surface. Properties of Porous Media. Bentheim sandstone samples were used because of the homogeneity, uniform porosity (22−24%), and high permeability (1.1D) of this medium. This sandstone has a grain density of 2.65 g/cm3, and the mineralogy typically showed 95−99% quartz content with traces of the clay mineral kaolinite. Table 1 shows some of the rock specific parameters used for the numerical part of the study, and Table 2 lists the sample geometry and the saturations used in the experimental and numerical part. Additionally, Evolving Porous Media (EPM) option number 332 was invoked to account for changes in relative permeability as a function of reduced effective porosity due to hydrate formation. Experimental Procedure. Each sample was positioned between two circular open volume polyoxymethylene (POM) spacers to prevent water transport into the lines and subsequent hydrate plugging (see Figure II, Supporting Information, for picture of sample and core assembly). The core was mounted between two polyetheretherketone (PEEK) end-pieces with Teflon shrink tubing, and was positioned inside the composite fiberglass core holder. The sample was evacuated and partially saturated (44−50%) with 0.1 wt % NaCl brine solution to prevent clay swelling. The sample was assembled in the pressure vessel, in which we incrementally increased the pore- and overburden pressure to 8.38 and 10.44 MPa, respectively. The CH4 pump was then set to operate in constant pressure mode at 8.38 MPa, pressurizing the sample from both ends, while the system was cooled to 4 °C. CH4 in hydrate has higher density relative to the gas phase at the experimental conditions. Changes in hydrate saturation could therefore be quantified by monitoring changes in pump volume, where initiation of hydrate formation was indicated by sudden changes in pump volume and temperature, as well as by loss of the MRI signal. Pressure was reduced to 4.2 MPa using constant volumetric rate, and maintained above the equilibrium pressure (Peq = 4.02 MPa at 4 °C) at 4.2 MPa when no further variations in the MRI signal or the
Table 1. Rock Specific Parameters Used during the Numerical Study grain density [g/mL], ρgrain porosity [fraction], φ absolute permeability [m2], K dry thermal conductivity [W/m/K], kdry wet thermal conductivity [W/m/K], kwet composite thermal conductivity model capillary pressure model39
2.650 0.23 1.10 × 10−12 0.5 3.1
k Θ = kdry + ( SA +
SH )(k wet − kdry ) + ϕSIλI
Pc = − P0[(S*)−1/ λ − 1]1 − λ
S* = van Genuchten soil param., λ irreducible aqueous saturation, SirA air entry pressure, kPa max. pressure, Pmax max. aqueous saturation, SmxA relative permeability model (Mmodified version of Stone40)
(SA − SirA) (SmxA − SirA)
0.45 0.145 12.5 1.00 × 1006 1
⎧ ⎧⎡ (S − S ) ⎤n⎫⎫ ⎪ ⎪ irA ⎨⎢ A k rA = max⎨ 0, min ⎥ ⎬⎬ ⎪ ⎪ ⎩⎣ (1 − SirA) ⎦ ⎭⎭ ⎩ ⎪
⎪
⎪
⎪
⎧ ⎧⎡ (S − S ) ⎤nG ⎫⎫ ⎪ ⎪ irG k rG = max⎨ 0, min⎨⎢ G ⎥ ⎬⎬ ⎪ ⎪ ⎩⎣ 1 − SirA ⎦ ⎭⎭ ⎩ SirA SirG aqueous exponent, n gas exponent, nG
⎪
⎪
⎪
⎪
0.15 0.02 3.5 3.2
pump volume were detected. No intensity or volume changes were observed as the system was given time to equilibrate. The right core boundary was isolated from the pumps by closing a valve, while the left core boundary was kept open. Flow was therefore only permitted across the left core boundary. We then reduced the pressure to 3.96 MPa, which we maintained constant during dissociation. Only minor dissociation was observed before the system stabilized, which was attributed to water freshening or dissociation induced cooling. The inlet pressure was therefore further reduced to 3.89 MPa to promote further dissociation. Again, hydrate dissociation was retarded prior to complete dissociation, and the pressure was further reduced to 3.82 MPa to fully dissociate all gas hydrate in the sample. Tests 1 and 2 5575
dx.doi.org/10.1021/ef500255y | Energy Fuels 2014, 28, 5573−5586
Energy & Fuels
Article
Table 2. Sample Properties and Saturation Data Used for the Numerical Study Are Based on Data from Experiments
test 1 test 2 test 3
core length [cm]
core diam. [cm]
SA,initial
SA
SH
SG
15.0 15.0 14.2
3.81 3.81 5.06
0.45 0.50 0.44
0.045 0.050 2.6 × 10−4a
0.51 0.57 0.56
0.44 0.38 0.34
Table 3. Thermal Properties Assigned to Each Material in the Numerical Study thermal conductivity [W/m·K] specific heat [J/kg·K] Bentheim Fluorinert Teflon shrink tubing POM spacer PEEK endpiece CH4 SS316 Buna-N sleeve
a
Model was initialized based on a T+H equilibration process, where water was crystallized until the formation brine concentration reached 14 wt % NaCl. Experimental PVT data suggest less conversion, where capillary forces will inhibit conversion of water in smaller pores and along mineral grains. However, numerical data was emphasized in the simulation of test 3.
3.1a 0.065 0.25
1000 1100 1200
0.23 0.25 0.0346b 16.2 0.25
1465 1320 2260 502 1350
a
Values used in this study are for Berea sandstone, which has similar properties as Bentheim. bThermal conductivity at 3.95 °C and 3.82 MPa.
therefore have hydrate dissociation at three consecutive pressure steps (3.96, 3.89, and 3.82 MPa). A similar procedure was repeated for test 3, where pressure was reduced from 8.3 to 4.2 MPa after hydrate formation. After being equilibrated, pressure was reduced in a single step (at constant volumetric rate) to 3.2 MPa in order to investigate the dissociation behavior under conditions of a larger depressurization-induced driving force. Pressure was maintained constant at 3.2 MPa during hydrate dissociation for test 3.
order to accurately describe conductive heat exchange with the constant temperature boundary. Uniform discretization was used to represent the porous media, while the element size increased progressively with the distance from the vertical sides of the core. Heterogeneities may have a significant impact on the production rate,41 so the variable saturation distribution data from the MRI observations was accurately represented (to the extent possible) in the initial numerical model. As indicated in Figure 5, the outer radial positions of the sample was less hydrate saturated, where gas was occupying the remaining pore space. Voronoi MESH. Irregular Voronoi tessellation was required to fully take advantage of MRI data and maintain well-defined boundaries for accurate description of heat transfer. Within the hydrate core, the Voronoi mesh used the discretization provided from the MRI data as tessellation generator points. Consequently, the interior of the core was primarily composed of a Cartesian (i.e., rectilinear) mesh. The outer edges of the core, and the surrounding experimental apparatus, used a series of cylindrical rings of points as tessellation generator points. The mesh was thus composed of an inner Cartesian region, with properties populated by the MRI data, and an outer pseudocylindrical region comprising the experimental apparatus. A tetrahedral decomposition of one quadrant of the Voronoi mesh is depicted in Figure 6. This approach resulted in an overall reduced number of elements and connections, as compared to a fully discretized Cartesian mesh for the entire domain. The Voronoi mesh ultimately consisted of 137 054 elements with 537 331 connections, and therefore required significant computation time to solve the problems.
■
NUMERICAL DESCRIPTION Numerical Model. The physical system was described as a 2D Cartesian system and as a complex 3D Voronoi grid. The 2D grid was used for the majority of this study because of the significant computational and executional requirements of the 3D model, which involved the solution of 600.000 coupled equations of mass, heat, and flow at each time-step. Figure 4 illustrates the numerical representation of the physical system shown in Figure 3. The thin layer of Fluorinert surrounding the system represented a no-flow (Neuman) boundary, and the temperature measured by the thermocouple was used as a timevariable temperature boundary. The specific heat and thermal conductivity values of each material involved in this study are listed in Table 3. The left open spacer volume (gas filled, depicted by the white space within the POM spacer) was connected to a pump that maintained constant pressure, which was represented by a constant pressure element (Dirichlet-type boundary condition) in the numerical model. The position of this sink is illustrated by VALVE in Figure 4. Discretization of the Numerical Models. The saturation distribution data from the MRI observations were used to design the numerical mesh (see Figure 5). The 2D grid comprised 8765 elements, of which 2184 elements represented the core sample. The additional elements were necessary in
Figure 4. Illustration of the numerical system based on the cross-sectional view of the core holder in Figure 3. A no-flow (Neuman) boundary confined a thin layer of Fluorinert, and the thermocouple positioned in the confining fluid provided data for a time-variable temperature boundary. Specific heat and thermal conductivity of the different materials are provided in Table 3. 5576
dx.doi.org/10.1021/ef500255y | Energy Fuels 2014, 28, 5573−5586
Energy & Fuels
Article
Figure 5. Gas hydrate saturation distribution based on MRI. The figure clearly shows heterogeneities in the initial distribution, where higher concentrations of gas hydrate are present in the vicinity of the end-pieces, while the gas saturation is higher at outer radial positions close to the left and right end piece. These gas-saturated parts of the sample may adversely impact heat transfer within the sample.
Figure 6. Quadrant-view of a tetrahedral decomposition of the Voronoi mesh. The Voronoi mesh used a Cartesian discretization for the region representing the hydrate core, which transitioned into a pseudocylindrical region outside of the core. The cylindrical region precisely captures the locations of the various components of the experimental apparatus.
■
EXPERIMENTAL RESULTS
the constant pressure operating pumps (8.38 MPa). Good agreement between MRI intensity and pump volume changes during hydrate growth is demonstrated in Figure 7 for test 1. The MRI intensity is proportional to the amount of water in the system and was therefore inverted and normalized for easier comparison with CH4 consumption during growth. Small
Hydrate formation. Hydrate formation was detected as sudden increase in CH4 delivery from the pumps, as the CH4 density increases in the hydrate crystal relative to the gas phase at the experimental conditions. Changes in hydrate saturation were therefore quantified based on changes in pump volume for 5577
dx.doi.org/10.1021/ef500255y | Energy Fuels 2014, 28, 5573−5586
Energy & Fuels
Article
fugacity difference f v − feq). These two experiments (tests 1 and 2) also required three consecutive pressure steps (3.96, 3.89, and 3.82 MPa) for full dissociation, where the pressure was reduced when production ceased at a given pressure. Limited dissociation was observed at 3.96 MPa (5−9% hydrate mass) for tests 1 and 2. Constant pressure/temperature (up to ∼0.2 °C variation in the confining fluid) conditions were maintained for up to 80 h with only minor variations in MRI intensity and CH4 production, and the pressure was therefore reduced to 3.89 MPa to trigger further dissociation. The dissociation rate at 3.89 MPa was variable, with a deflection point after 90 h for both test 1 and test 2. This corresponds to hydrate reformation and may be a result of dissociation-related water freshening (thus causing a shift in the P−T equilibrium) and locally reduced temperatures. The dissociation trend was consistent for the two experiments; thus, we suspect that the observed behavior is related to low driving force and sudden changes in core temperature. In both experiments, the majority of the gas hydrate dissociated at 3.82 MPa. Dissociation did not occur spontaneously as the pump pressure was reduced to 3.82 MPa in test 2. The boundary temperature (Tb = ∼3.7−3.8 °C) was only slightly higher than the equilibrium temperature (Teq = 3.61 °C for 0.6 wt % NaCl brine concentration); however, it is reasonable to assume that the sample temperature was slightly lower. Also, we have to keep in mind that previous research42−45 have reported cases where gas hydrates are preserved for extended periods above their hydration temperature. Still, retention in hydrate dissociation may simply be related to equilibrium conditions induced by a temperature drop (Figure IV, Supporting Information). Temperature Effects during Dissociation. Hydrate dissociation is a strongly endothermic process. This selfpreserving mechanism will result in regained stability following a dissociation event unless heat is provided through advection or conduction. The dissociation rate in our experiments was reduced by dissociation-induced temperature fluctuations, which shifted the boundary temperature with a subsequent shift in the hydration pressure. The driving force may also be expressed in terms of subcooling, where
Figure 7. Variations in MRI intensity and CH4 consumption corresponded well during hydrate formation. Confining temperature was maintained constant at 4.0 °C toward the end of hydrate formation. The hydrate formation process mainly occurred during the initial 20 h, but the temperature/pressure conditions were maintained constant for additional 140 h to allow the system to equilibrate. Only minor changes were observed during this period.
temperature variations were observed during hydrate growth due to the exothermic reaction, but the temperature was maintained at 4 °C toward the end of hydrate formation and at the initiation of the depressurization process. Mass Balance during Hydrate Dissociation. Three different depressurization experiments were conducted, where two experiments (tests 1 and 2) were assisted by in situ MR imaging, and will therefore be emphasized in the discussion. Dissociation was quantified based on changes in pump volume, which were consistent with changes in MRI intensity (see Figure III, Supporting Information). Hydrate dissociation results from all three experiments are included in Figure 8. Significant differences were observed in terms of dissociation time, where complete dissociation required less than 30 h for test 3 (induced by Pb = 3.2 MPa). In comparison, more than 280 h were required when dissociation occurred at a higher pressure (tests 1 and 2), which minimizes the dissociation driving force (determined by the
Subcooling = Tequilibrium − Tsample
Positive subcooling indicates a system within the hydrate stable region, while negative subcooling will trigger hydrate dissociation. Tests 1 and 2 were initially characterized by low negative subcooling (T1 = −0.11 K; T2 = −0.06 K (the degree of subcooling assumes temperature equilibrium between the measured confining fluid temperature and sample temperature)), while test 3 had subcooling of −2.17 K. We therefore expected a different depressurization response for the three experiments. Gas hydrate dissociation is compared to subcooling in Figure 9. Equilibrium temperatures were estimated from the statistical thermodynamics code CSMGem,46 while a moving average function was used to better visualize trends in temperature fluctuations (due to the 0.1 °C resolution). Experimentally observed inflection points in CH4 production rate corresponded well with variations in subcooling. The degree of subcooling increased (toward 0) as gas hydrates dissociated, induced by a combination of water freshening and temperature variations in the confining fluid. Finally, gas hydrate dissociation was retarded as the subcooling reached 0. The rapid temperature increase (reflected in the subcooling response)
Figure 8. Gas hydrate dissociation was quantified based on pump volume changes and was normalized to moles of CH4 consumed during hydrate growth. Gas hydrates were dissociated at three pressure steps, where dissociation trends at each of the pressure steps were more or less consistent for tests 1 and 2. Dissociation time was reduced by an order of magnitude at lower dissociation pressure (test 3). 5578
dx.doi.org/10.1021/ef500255y | Energy Fuels 2014, 28, 5573−5586
Energy & Fuels
Article
Figure 9. Experimental conditions in test 1 were maintained in the vicinity of the three-phase line (Lw−V−H), where small temperature variations affected the degree of subcooling. Observed inflection points in dissociation trends corresponded well with changes in subcooling, where regained stability (or hydrate reformation) was observed for positive subcooling.
Figure 10. 1D saturation profiles along the core during dissociation in test 1. Two profiles are included at each pressure step for comparison of initial (light colors) and final (darker colors) saturation state. The profiles demonstrate how heat transport across the end-pieces enhances the local dissociation, illustrated by the increased Sw. The left and right intensity peaks are due to gas/water accumulation in the POM spacer volume (originally gas filled).
observed after 280 h may indicate complete hydrate dissociation in the vicinity of the thermocouple.24 However, similar temperature increase was not observed for test 2 toward the end of dissociation (see Figure IV, Supporting Information). Dissociation may be limited to the water/hydrate interface for low-permeability sediments. Limited pressure propagation into the core may have a significant impact on the dissociation performance, resulting in localized pressure and temperature conditions that can remain within the hydrate stable region. We recorded pressure drops of up to 124 kPa (1.24 bar) across the core length during the pressure reduction phase of our experiment, but recorded differential pressures were usually below 6.9 kPa (0.07 bar) during the dissociation operation; thus, they were not expected to be of major importance. Saturation from MRI during Dissociation. The mass balance data were complemented by the analysis of continuous MR imaging during dissociation, either as 1D profiles (test 1) or as full 3D images (test 2). Figure 10 indicates that the initial dissociation rate was more pronounced in the vicinity of the end-pieces (i.e., near the POM spacers) than in the middle section of the core in test 1. We believe heat transfer was the dominant mechanism, as only the left transaxial core face was connected to the pump during production, while intensity buildup was observed at both core faces (x = 0.3 cm and x = 15.3 cm). Figure 10 indicates that the average Sw at the core faces exceeded 0.3 at 3.89 MPa, while the middle core segment Sw was less than 0.2. Hydrate dissociation also occurred along the cylindrical surface at the core periphery where heat was transported across the Teflon shrink tubing. Residual water persisted in the right midsection of the test 2 sample, despite significant subcooling during formation (Figure 11a and x = 8−13 cm in Figure V, Supporting Information). The gas hydrate stability regime may have been affected by locally higher NaCl concentrations, or by water under-saturated in CH4. Diffusivity of CH4 in hydrate is a slow process, and water occluded by thicker hydrate layers may have remained under-saturated in CH4. Dissociation at Pb = 3.96 MPa was mainly observed adjacent to the residual water, where mobilization of nonfrozen water may have altered the local hydrate stability (Figure V, Supporting Information). This may
be attributed to nonequilibrium conditions imposed by the excess water. Dissociation at the Pb = 3.89 MPa pressure was focused in the proximity of the end-pieces, while dissociation progressed more uniformly when the pressure was Pb = 3.82 MPa (see Figure VI and VIII, Supporting Information). Hydrate dissociation and reformation was observed to occur as parallel events at Pb = 3.89 MPa, where water crystallized in the middle section of the sample while dissociation was observed adjacent to the end pieces. This was reflected in both the material balance and the MRI observations (Figures III and VII, Supporting Information) and was assumed to be due to a temperature gradient as no differential pressure was observed across the sample. The variable dissociation rates corresponded well with temperature changes in the confining fluid and demonstrate the importance of heat transport during dissociation. Production of Associated Water. Some production of associated water was observed during the dissociation process in test 2 (Figure IX, Supporting Information). This reduced the average MRI intensity measured across the sample, with a subsequent slight mismatch between the MRI intensity and CH4 production (Figure III, Supporting Information). A significant intensity increase was observed in the left (front) spacer volume after 96 h, which was attributed to water production. Water presence persisted during the dissociation process, as dissociation induced gas flux from the sample will minimize the potential for spontaneous imbibition. Water presence in the spacer volume was reduced as gas hydrates were reforming in the sample (after 167 h in Figure IX), where water imbibed back into the sample due to capillary forces. Water imbibition from the spacer volume was also observed after complete hydrate dissociation (after 352 h, Figure IX, Supporting Information).
■
NUMERICAL RESULTS Performance of the Equilibrium and the Kinetic Reaction Models. As a matter of simplicity and convenience, the equilibrium model is generally preferred because it involves one less equation per element. Kowalsky and Moridis33 5579
dx.doi.org/10.1021/ef500255y | Energy Fuels 2014, 28, 5573−5586
Energy & Fuels
Article
Figure 11. Water saturation from MRI during dissociation. Hydrogen in gas hydrate is not visible in the MRI, and signal therefore indicates presence of free water. Dissociation was focused close to the end-pieces and at the outer radii, where heat transport through conduction is a dominant mechanism. (a) Initial saturation state at 3.96 MPa, where dissociation occurred in the vicinity of the excess water. (b) Water mobilizing toward the left producer at 3.89 MPa. (c) Progression of dissociation close to the end-pieces with parallel hydrate growth in the middle section of the core. (d) More or less uniform saturation at the end of dissociation. The cut-out is for illustration purposes only to better see the saturation distribution within the sample.
accompanied by limited time-step size and longer simulations. Periodic drops were observed in the numerical prediction of the volumetric production in test 3 (see Figure XI, Supporting Information). These drops are consistent with observations from Kowalsky and Moridis.33 In comparison, the kinetic reaction model was less sensitive to temperature variations at the boundary, and required fewer time steps and shorter execution times to solve the problem.
demonstrated how the equilibrium and kinetic reaction model responded similarly during long-term (multiyear) dissociation, and reported minor differences at the reservoir scale. However, these differences were more pronounced on smaller scale domains (1.5 m) and for short-time processes (such as the dissociation of hydrates in core samples). This clearly showed the possible importance of kinetics in short-term processes (conforming to expectations) and indicated the possibility that the choice of reaction model depends on the physical problem size, especially when the small size leads to fast dissociation. A preliminary set of simulations of simplified representations of the laboratory investigations in test 1 involved a coarse grid (830 subdomains, 1753 connections) and resulted in similar predictions from the equilibrium and kinetic reaction model (Figure X, Supporting Information). The equilibrium model is sensitive to the P−T conditions during dissociation, small variations of which may induce abrupt changes and fluctuations in the production rate that are caused by treating dissociation as an equilibrium process. As hydrates dissociate at a given location, the resulting temperature drop and pressure increase or stability (as gas is released) adversely affect further dissociation of the hydrate in the vicinity, resulting in a temporary drop in production. Production increases again after sufficient heat flowing from the boundaries becomes available, resulting in this oscillatory behavior, which is often
■
NUMERICAL REPRESENTATION OF GAS HYDRATE DISSOCIATION Effect of Different Kinetic Combinations. We used 17 different combinations of kinetic parameter values in the preliminary part of this study in order to minimize the deviations between the numerical predictions and the laboratory observations. The activation energy (Ea) is the minimum energy that must be available in the system for a reaction to occur. As such, higher activation energy will adversely impact the predicted hydrate dissociation rate, while increasing intrinsic reaction constant and area adjustment factor will increase the hydrate dissociation rate. Different sets of parameter combinations may yield solutions similar to the experimental observations, and we therefore wanted to identify which combinations provided best predictions. Some combi5580
dx.doi.org/10.1021/ef500255y | Energy Fuels 2014, 28, 5573−5586
Energy & Fuels
Article
conclusion that the shrink tubing was properly fitted around the sample. Dynamic Kinetic Parameters. Dynamic (that is, varying during the simulation, as opposed to being constant) kinetic parameters may be necessary to account for variable hydrate morphology and interaction during growth/decomposition. Simulation 4−4 is an example of curve fitting, where kinetic parameters were adjusted twice during the simulation (4−4a, b, and c in Table 4. See also Figure XII, Supporting Information). This resulted in excellent agreement between empirical and numerical data. Effect of Heterogeneous Fluid Distribution. Heterogeneity in the fluid saturation distribution affects both the relative permeability and the heat transfer, and was initially presumed to play a significant role in the dissociation behavior. We therefore compared results in test 1 using a uniformly saturated model (HOM: SH = 0.51 in each element, as listed in Table 2) with a nonuniformly saturated model (4−11: see saturation distribution in Figure XIII, Supporting Information), both models employing the kinetic parameters of model 4−11. Laboratory observations suggested that 9.4% of the hydrate mass initially in place dissociated at Pb = 3.96 MPa (test 1). In comparison, numerical predictions utilizing the 4−11 parameter combination suggested 11.2% of the hydrate mass dissociating at this Pb, a prediction that was similar for both the uniformly and nonuniformly saturated sample (purple and green dashed line in Figure 12). Saturation differences became increasingly more pronounced during dissociation at Pb = 3.89 MPa. Figure 13 demonstrates a
nations used in this study are listed in Table 4 and illustrated in Figure 12. Combinations 4−0 (ΔEa = 8.1 × 104 kJ/mol K0 = Table 4. Kinetic Parameters Applied during Simulation of Empirical Data from Test 1 simulation file 4−0 4−4a 4−4b 4−4c 4−8 4−11 26by7
Ea [J/mol] 8.10 8.50 8.50 8.50 8.30 8.50 8.10
× × × × × × ×
104 104 104 104 104 104 104
intrinsic rate [mol/m2·Pa· s] 3.60 3.60 3.60 3.60 5.00 7.00 3.60
× × × × × × ×
104 105 104 104 104 104 104
area adjustment factor [dimensionless] 6.50 3.65 3.65 6.50 3.65 3.65 3.65
× × × × × × ×
10−2 10−1 10−1 10−2 10−1 10−1 10−1
gas hydrate thermal conductivity [W/m/K] 4.50 4.50 4.50 4.50 5.70 5.70 4.50
× × × × × × ×
10−1 10−1 10−1 10−1 10−1 10−1 10−1
Figure 12. Comparison of experimental and numerical dissociation results at 3.96 MPa for test 1. A range of kinetic parameters were used for increased prediction accuracy. Best match was achieved through simulation 4−0 and 4−11, where the kinetic values are listed in Table 4.
3.6 × 104 mol/m2··Pa·s, FA = 0.065) and 4−11 (ΔEa = 8.5 × 104 kJ/mol, K0 = 7.0 × 104 mol/m2·Pa·s, FA = 0.365) were identified as the ones with the optimal parameter values, as they resulted in the best agreement between predictions and laboratory observations. These parameter values were employed for the remainder of this study. Comparison of 2D and 3D Predictions. The computational difficulty is dependent upon the spatial discretization of the domain. Reducing a 3D problem by one dimension will significantly reduce the computational cost; however, the model performance should be validated through comparison. A 3D Voronoi model was used to verify the performance of the 2D model predictions, where two different initialization strategies were used: (1) the 3D Voronoi model was initialized with fluid saturation according to MRI data (V_HET), or (2) the 3D Voronoi model was initialized with homogeneous fluid saturation, but the model assumed the existence of a thin gas layer between the sample and Teflon shrink tubing (V_Gas). The kinetic parameter combination 4−0 was used for both 3D scenarios, and the predictions are visualized in Figure 12. Simulation 4−0 (2D) and V_HET (3D) utilized identical input parameters and exhibited similar trends, thus validating and justifying the 2D simplification. V_Gas (3D) did not reproduce the observations in a satisfactory manner, leading to the
Figure 13. Deviations between nonuniformly and uniformly saturated samples in terms of dissociation rates were less dominant than first assumed. Images compare water saturation at 3.89 MPa after 81.3 h, where predictions deviated by 4% in saturation units. Dissociation was more pronounced for elements in contact with the POM material (corner elements), as heat conductivity in POM is 1 order of magnitude higher relative to CH4.
4% difference in average SH for the two predictions at t = 81.3 h. The fluid distribution in the uniformity-based model is relatively homogeneous, while the spatial saturation variations in the heterogeneous model still resemble the initial experimental conditions depicted in Figure 5. This indicates that kinetic dissociation progresses in a relatively uniform manner, which is consistent with core-scale observations,33 but also emphasizes that the dissociation process is affected by saturation heterogeneity, which may be related to saturation dependent heat and fluid transport. 5581
dx.doi.org/10.1021/ef500255y | Energy Fuels 2014, 28, 5573−5586
Energy & Fuels
Article
relative to test 1 (5% vs 9%) at Pb = 3.96 MPa. This was attributed to lower temperature boundary in test 2, thus reducing the dissociation driving force. Hydrate dissociation and subsequent reformation was observed at Pb = 3.89 MPa in both empirical and numerical trends. This process was linked to boundary temperature variations, where ∼0.10 °C variation shifted the process from hydrate dissociation to reformation. The temperature impact in terms of dissociation rate was consistently more pronounced in the 4−11 data set due to higher intrinsic rate. In Figure 15, the predicted hydrate dissociation (in mass fraction) was within 3% of experimental observations for Pb = 3.96 MPa. The deviation was more pronounced for Pb = 3.89, in which case the 4−11 data set overestimated mass fraction of dissociated hydrate by 16%. In comparison, the final dissociation estimate by data set 4−0 was within 0.4% of the experimentally observed mass fraction. Both numerical predictions exhibited an immediate response to the temperature drop after 77 h (Figure 15), which was accompanied by hydrate reformation in the middle core segment. The laboratory data of dissociation at Pb = 3.82 MPa were characterized by three temperature plateaus, which dominated the dissociation rate. This was not as well reflected in the numerical simulation results, but predictions were consistently within the same range as the laboratory data. This may be attributed to inconsistency between the fluid and medium thermal equilibrium assumption. In the model, we assume thermal equilibrium in every element; that is, we assume that the fluid and matrix elements have the same temperature. This assumption will only be true after a long time, when such thermal equilibrium has been established. Otherwise, the measured temperature is likely lower than the numerically predicted one, with subsequent overprediction of hydrate dissociation. Spatial Saturation Variations during Predicted Dissociation. Figure 10 demonstrated through 1D saturation profiles how experimental dissociation was slightly more abundant in vicinity of the transaxial core faces, probably resulting from favorable heat transport conditions. This was also observed in numerical predictions, as illustrated in Figure 16. Water was more abundant in vicinity of the transaxial core faces (horizontal position 0 and 15 cm), where average element temperature was higher due to favorable heat transport conditions (i.e., presence of POM plastic). Predicted temperatures were generally higher at the left transaxial core face (horizontal position = 0 in Figure 16) with minimum temperatures in the middle core segment. Initial dissociation in Figure 16 (t = 3.45 h) was to some extent limited by pressure propagation into the sample. Abundant dissociation in vicinity of the constant pressure sink (at horizontal position = 0) resulted in reduced temperature and locally elevated pressures due to hydrate dissociation. This may have provided favorable conditions for gas hydrates to reform and may explain the spatial Sw trends observed after 10 h. Saturation differences depicted after 110 h shows the most commonly predicted saturation trend, where Sw was more abundant in elements with favorable heat transport. The thermal conductivity of POM is 1 order of magnitude higher than that of CH4 (Table 3), and thus, the POM material will dominate heat conduction across the core face. Gas hydrates were slightly less abundant at radially distant positions, which was consistent with observations from MRI. Overall, heat
A uniform initial saturation distribution was computationally more desirable, as it required shorter execution times. The MRI visualization data in test 1 were obviously more relevant as they represented the actual system. The initial saturation distributions in test 2 increased markedly the computational difficulty of the problem because significant excess water in the middle segment of the core and an imperfect space discretization resulted in significant saturation variations between adjacent gridblocks, thus complicating convergence and reducing the time step size. Convergence was not achieved unless a more uniformly saturation distribution was applied (due to the nonequilibrium condition of the problem). Test 2 was therefore initialized with uniform saturation distribution in combination with kinetic parameters from the 4−0 and 4−11 parameter data sets. Kinetic Predictions of Dissociation: Test 1. In Figure 14, we compare laboratory data and numerical predictions
Figure 14. Comparison of experimental and numerical dissociation results for test 1 at three different pressure steps. The general agreement was good, but the numerical data were less sensitive to temperature variations in the confining fluid. Variations in dissociation rate were observed, but neither 2D nor 3D modeling was capable of reproducing empirical reformation for test 1.
associated with the depressurization experiment (test 1) induced by Pb = 3.96, 3.89, and 3.82 MPa. As expected, the dissociation parameter data set with the higher intrinsic rate (4−11) overestimated hydrate dissociation, while predictions with the 4−0 parameter data set were within the same range as the laboratory observations. The dissociation response is also affected by the activation energy, but the intrinsic rate was dominant as it was almost twice as large for combination 4−11 relative to combination 4−0. None of the numerical predictions managed to reproduce the inflection point observed toward the end of the depressurization experiment with Pb = 3.89 MPa, which we attributed to hydrate reformation. This may indicate that the initialization of the numerical model was imperfect, or that temperature gradients in the samples could not be accurately described because of insufficient number of thermocouples (one only), or that there was variability (spatial and/or temporal) in the parameters that had not been captured in the model. Kinetic Predictions of Dissociation: Test 2. Predicted and observed dissociation in test 2 at all pressure steps are summarized in Figure 15. Both models were initialized with homogeneous fluid distribution according to Table 2, while kinetic combinations 4−0 and 4−11 were used. A smaller mass fraction of gas hydrate dissociated in the test 2 experiment 5582
dx.doi.org/10.1021/ef500255y | Energy Fuels 2014, 28, 5573−5586
Energy & Fuels
Article
Figure 15. Comparison of empirical and numerical hydrate dissociation for test 2 at 3.96, 3.89, and 3.82 MPa. Small variations in boundary temperature affected experimentally observed and numerically predicted hydrate dissociation, as illustrated at 3.96 MPa. Numerical deviations were more pronounced at 3.89 MPa but were overall consistent with empirical trends.
therefore sensitive to small deviations in ion concentration, as illustrated in Figure 17, where numerical predictions utilizing salt concentrations ranging between 0.008 and 0.010 (by mass fraction) have been included. Dissociation induced salt dilution will not trigger hydrate reformation independently, but the system becomes more sensitive to temperature variations as the dissociation driving force decreases (due to water-freshening). Temperature Sensitivity during Hydrate Dissociation. Gas hydrate stability and dissociation is sensitive to temperature variations induced by the endothermic reaction occurring within the sample or from outside stimuli (such as variations in
transport appeared to be the dominant mechanism for experimentally and numerically observed dissociation. Sensitivity to Variations in Initial Salt Concentration. Na+ and Cl− interact with water molecules and affect the gas hydrate stability. The initial NaCl concentration for test 1 was 0.1 wt % (before hydrates were formed). The salt concentration increases during hydrate formation because salt cannot become a part of the hydrate crystal and is thus expelled into the remaining water, and was estimated to be 0.9 wt % (0.009 by mass fraction) at the completion of growth for test 1. The dissociation test was conducted at low driving force and was 5583
dx.doi.org/10.1021/ef500255y | Energy Fuels 2014, 28, 5573−5586
Energy & Fuels
Article
Figure 18. Gas hydrate dissociation rate was minimized as dissociation induced water-freshening shifted the sample conditions toward the thermodynamic equilibrium. Immediate instability and reformation occurred in response of the variation in boundary temperature after 4.6 days, where temperature was reduced by either −0.01 °C, −0.05 °C, or −0.10 °C. Figure 16. Initial decomposition was more abundant at horizontal position = 0 (illustrated by abundance of water), where heat and fluids were advancing toward the constant pressure sink. The last image (t = 110 h) demonstrates how dissociation is facilitated in corner elements due to favorable heat conductivity, which corroborates empirical observations. Color bars have different scale to highlight saturation differences at each time-step.
Equilibrium and Kinetic Predictions at Higher Dissociation Driving Force. Sequential depressurization at higher pressures required more than 280 h for complete dissociation. In comparison, less than 30 h were required when dissociation occurred at Pb = 3.2 MPa in test 3. Empirical and numerical predictions at 3.2 MPa, employing both the equilibrium and kinetic reaction model, are compared in Figure 19. The kinetic variables of the 4−0 data set resulted in
Figure 17. Comparison of numerically predicted dissociation at 3.96 MPa for test 1 using different salt concentrations (0.008 to 0.01 by mass fraction). Small deviations in NaCl concentration impacts the dissociation rate. Correct initialization of system properties is therefore essential for correct prediction.
Figure 19. Comparison of empirical and numerical dissociation at 3.2 MPa. Kinetic combination 4−0 resulted in almost identical behavior as experimental observations, while the equilibrium model significantly overestimated dissociation. The empirical dissociation rate was sensitive to minor temperature variations, despite significant driving force. This confirms that heat transfer is a controlling mechanism during dissociation. T1 = sample temperature, T2 = boundary (confining) temperature.
the boundary temperature). Temperature sensitivity was investigated in a 2D model, where hydrate dissociation occurred at Pb = 3.96 MPa, with a constant boundary temperature of 3.95 °C. Dissociation occurred over several days due to low driving force and was eventually retarded due to water freshening. An immediate reformation response was observed in Figure 18 as the boundary temperature was shifted by either −0.01 °C, −0.05 °C, or −0.10 °C after 4.58 days, thus impacting the thermodynamic equilibrium. The extent of hydrate reformation was more pronounced at lower temperatures and therefore emphasizes the importance of using an accurate time-variable temperature boundary when small temperature variations are observed in the experiments.
excellent agreement with the laboratory measurements, while the equilibrium reaction model underestimated the dissociation time by a factor of 6 (30 h versus 5 h). Temperature sensitivity was maintained despite significant driving force, demonstrated by the kink in dissociation rate after 9.5 h. The choice of reaction model should therefore be reflected by the physical size of the problem, where small-scale rapid dissociation may require a kinetic reaction model due to kinetic limitations. In this study, deviations in the mass accumulation term were more 5584
dx.doi.org/10.1021/ef500255y | Energy Fuels 2014, 28, 5573−5586
Energy & Fuels
Article
Author Contributions
pronounced at higher driving force, in conformity with our expectations for short-term dissociation.
⊥
The manuscript was written through contributions of all authors. All authors have given approval to the final version of the manuscript. K.A.B., C.M.F., G.J.M., and A.G. contributed equally.
■
CONCLUSIONS The main scope of this study was to determine the accuracy of numerical predictions relative to empirical data. Three dissociation experiments on partially hydrate saturated samples confirmed that heat transport was one of the main limitations during dissociation, where the hydrate dissociation rate was shifted by minor boundary temperature variation. This was observed irrespective of initial driving force, but the impact was more pronounced for the two experiments with low negative subcooling, where hydrate reformation was triggered by a combination of reduced temperatures and water freshening. Gas hydrates are complex by nature, and an MRI therefore assisted in monitoring phase transitions within the porous media. MRI saturation data corresponded well with mass balance data from pump logs, where dissociation appeared more abundant in vicinity of the end-pieces and at radial distant positions. This was later predicted numerically. The MRI revealed core segments with noncrystallized water after gas hydrates had formed. The reason for this is not known, but it may relate to hydrate blocking, where the water is occluded, or the water may have been under-saturated in CH4. Heterogeneous behaviors, both during hydrate formation and dissociation, clearly demonstrates the advantage of monitoring saturations and phase transitions in situ. Empirical data trends were overall well predicted numerically for all three dissociation experiments, where small changes in the confining fluid temperature were accounted for by implementing a time-variable boundary temperature. Inflection points in numerical predictions and experimental observations were related to varying boundary temperature. Consistency between predictions and observations instill confidence in the numerical code T+H and demonstrate the necessity of sufficient heat transport during hydrate dissociation. We have found that kinetic reaction modeling was necessary on core scale simulations, both in terms of accuracy and computational time. While equilibrium modeling is generally preferred, reaction kinetics may be needed to describe complex gas hydrate behavior within porous media, where fluid activity, capillary forces, contact area, mineral interaction, and other physical mechanisms impact the dissociation rate and results in nonequilibrium conditions. Deviations between the equilibrium and kinetic reaction model were more pronounced at higher driving force, which emphasizes the need for a kinetic reaction model for short-term dissociation. Overall good agreement was achieved in this study by utilizing kinetic parameters ΔEa = 81 kJ/mol, K0 = 3.6 × 104 mol/m2·Pa·s, and FA = 0.065.
■
Notes
The authors declare no competing financial interest.
■
ACKNOWLEDGMENTS Thanks to Dr. Jarle Husebø (from Statoil), Dr. Geir Ersland (from the University of Bergen), Jim Stevens and Dr. James Howard (from ConocoPhillips) for their contribution and experimental expertise. Thanks to Christian Hågenvik, Hans Berge, Truls Hamre Håheim, and Reza Hossainpour (University of Bergen) for their contribution on test 3. One of the authors is indebted to the Research Council of Norway, Statoil, and ConocoPhillips for funding.
■ ■
ASSOCIATED CONTENT
Figures I−XIII: Additional graphs. This material is available free of charge via the Internet at http://pubs.acs.org.
AUTHOR INFORMATION
Corresponding Author
*E-mail:
[email protected]. Phone: +47 45473748. Present Address ∥
REFERENCES
(1) Circone, S.; Kirby, S.; Stern, L. A. Direct Measurement of Methane Hydrate Composition along the Hydrate Equilibrium Boundary. J. Phys. Chem. B 2005, 109, 9468−9475. (2) Sloan, E. D.; Koh, C. A. Clathrate Hydrates of Natural Gases, 3rd ed. CRC Press: Boca Raton, 2008. (3) Moridis, G. J.; Collett, T. S.; Boswell, R.; Kurihara, M.; Reagan, M. T.; Koh, C.; Sloan, E. D. Toward Production From Gas Hydrates: Current Status, Assessment of Resources, and Simulation-Based Evaluation of Technology and Potential. SPE Reservoir Eval. Eng. 2009, 12, 745−771. (4) Wilder, J. W.; Moridis, G. J.; Wilson, S. J.; Kurihara, M.; White, M. D.; Masuda, Y.; Anderson, B. J.; Collett, T. S.; Hunter, R. B.; Narita, H.; Pooladi-Darvish, M.; Rose, K.; Boswell, R. an International Effort to Compare Gas Hydrate Reservoir Simulators. 6th International Conference on Gas Hydrates, Vancouver, BC, Canada; 2008. (5) Hong, H. N.; Pooladi-Darvish, M. A Numerical Study on Gas Production From Formations Containing Gas Hydrates. Petroleum Society of Canada: Canadian International Petroleum Conference, Calgary, AB, Canada: 2003. (6) Kurihara, M.; Kunihiro, F.; Ouchi, H.; Masuda, Y.; Narita, H. Investigation on Applicability of Methane Hydrate Production Methods to Reservoirs with Diverse Characteristics. 5th International Conference on Gas Hydrates, Trondheim, Norway; 2005. (7) Moridis, G. J. Numerical Studies of Gas Production From Methane Hydrates. SPE J. 2003, 8, 359−370. (8) Moridis, G. J.; Kowalsky, M. B., Gas Production from Unconfined Class 2 Hydrate Accumulations in the Oceanic Subsurface. In Economic Geology of Natural Gas Hydrates; Max, M. D., Ed; Kluwer Academic/Plenum Publishers: New York, 2005; pp 249−266. (9) Moridis, G. J.; Kowalsky, M. B.; Pruess, K. DepressurizationInduced Gas Production From Class 1 Hydrate Deposits. SPE Reservoir Eval. Eng. 2007, 10, 458−481. (10) Moridis, G. J.; Reagan, M. T. Strategies for Gas Production From Oceanic Class 3 Hydrate Accumulations. Offshore Technology Conference, Houston, TX; 2007. (11) Anderson, B. J.; Kurihara, M.; White, M. D.; Moridis, G. J.; Wilson, S. J.; Pooladi-Darvish, M.; Gaddipati, M.; Masuda, Y.; Collett, T. S.; Hunter, R. B.; Narita, H.; Rose, K.; Boswell, R. Regional LongTerm Production Modeling from a Single Well Test, Mount Elbert Gas Hydrate Stratigraphic Test Well, Alaska North Slope. Mar. Pet. Geol. 2011, 28, 493−501.
S Supporting Information *
■
ABBREVIATIONS T+H = TOUGH+HYDRATE MRI = magnetic resonance imaging POM = polyoxymethylene PEEK = polyetheretherketone
ConocoPhillips, Ekofiskvegen 35, 4056 Tananger, Norway 5585
dx.doi.org/10.1021/ef500255y | Energy Fuels 2014, 28, 5573−5586
Energy & Fuels
Article
Depressurization in Hydrate-Bearing Porous Medium. Energy Fuels 2012, 26, 1681−1694. (29) Sun, X.; Mohanty, K. K. Kinetic Simulation of Methane Hydrate Formation and Dissociation in Porous Media. Chem. Eng. Sci. 2006, 61, 3476−3495. (30) Tang, L.-G.; Li, X.-S.; Feng, Z.-P.; Li, G.; Fan, S.-S. Control Mechanisms for Gas Hydrate Production by Depressurization in Different Scale Hydrate Reservoirs. Energy Fuels 2006, 21, 227−233. (31) Yousif, M. H.; Abass, H. H.; Selim, M. S.; Sloan, E. D. Experimental and Theoretical Investigation of Methane−Gas− Hydrate Dissociation in Porous Media. SPE Reservoir Eng. 1991, 6, 69−76. (32) Moridis, G. J.; Kowalsky, M. B.; Preuss, K. TOUGH+HYDRATE v1.2 Users Manual: A Code for the Simulation of System Behavior in Hydrate-bearing Geologic Media; Lawrence Berkeley National Laboratory: Berkeley, CA, 2012. (33) Kowalsky, M. B.; Moridis, G. J. Comparison of Kinetic and Equilibrium Reaction Models in Simulating Gas Hydrate Behavior in Porous Media. Energy Convers. Manage. 2007, 48, 1850−1863. (34) Clarke, M.; Bishnoi, P. R. Determination of the Activation Energy and Intrinsic Rate Constant of Methane Gas Hydrate Decomposition. Can. J. Chem. Eng. 2001, 79, 143−147. (35) Kim, H. C.; Bishnoi, P. R.; Heidemann, R. A.; Rizvi, S. S. H. Kinetics of Methane Hydrate Decomposition. Chem. Eng. Sci. 1987, 42, 1645−1653. (36) Kvamme, B.; Graue, A.; Buanes, T.; Kuznetsova, T.; Ersland, G. Effects of Solid Surfaces on Hydrate Kinetics and Stability. Geological Society, London, Special Publications 2009, 319, 131−144. (37) HusebøJ.Monitoring Depressurization and CO2−CH4 Exchange Production Scenarios for Natural Gas Hydrates. Doctoral Thesis, University of Bergen, Bergen, Norway, 2008 (38) Husebø, J.; Graue, A.; Kvamme, B.; Stevens, J.; Howard, J. J.; Baldwin, B. A. Experimental Investigation of Methane Release from Hydrate Formation in Sandstone through Both Hydrate Dissociation and CO2 Sequestration. 6th International Conference on Gas Hydrates, Vancouver, BC, Canada; 2008. (39) van Genuchten, M. T. A Closed-form Equation for Predicting the Hydraulic Conductivity of Unsaturated Soils 1. Soil Sci. Soc. Am. 1980, 44, 892−898. (40) Stone, H. L. Probability Model for Estimating Three-Phase Relative Permeability. J. Pet. Technol. 1970, 22, 214−218. (41) Reagan, M. T.; Kowalsky, M. B.; Moridis, G. J.; Silpngarmlert, S. The Effect of Reservoir Heterogeneity on Gas Production From Hydrate Accumulations in the Permafrost. SPE Western Regional Meeting, Anaheim, CA; 2010. (42) Davidson, D. W.; Garg, S. K.; Gough, S. R.; Handa, Y. P.; Ratcliffe, C. I.; Ripmeester, J. A.; Tse, J. S.; Lawson, W. F. Laboratory Analysis of a Naturally Occurring Gas Hydrate from Sediment of the Gulf of Mexico. Geochim. Cosmochim. Acta 1986, 50, 619−623. (43) Stern, L. A.; Circone, S.; Kirby, S. H.; Durham, W. B. Temperature, Pressure, and Compositional Effects on Anomalous or “Self” Preservation of Gas Hydrates. Can. J. Phys. 2003, 81, 271−283. (44) Stern, L. A.; Kirby, S. H.; Durham, W. B. Peculiarities of Methane Clathrate Hydrate Formation and Solid-State Deformation, Including Possible Superheating of Water Ice. Science 1996, 273, 1843−1848. (45) Takeya, S.; Ebinuma, T.; Uchida, T.; Nagao, J.; Narita, H. SelfPreservation Effect and Dissociation Rates of CH4 Hydrate. J. Cryst. Growth 2002, 237−239 (Part 1), 379−382. (46) Ballard, A. L.; Sloan, E. D., Jr. The Next Generation of Hydrate Prediction: I. Hydrate Standard States and Incorporation of Spectroscopy. Fluid Phase Equilib. 2002, 194−197, 371−383.
(12) Kurihara, M.; Sato, A.; Ouchi, H.; NARITA, H.; Ebinuma, T.; SUZUKI, K.; Masuda, Y.; Saeki, T.; Yamamoto, K.; Fujii, T. Prediction of Production Test Performances in Eastern Nankai Trough Methane Hydrate Reservoirs Using 3D Reservoir Model. Offshore Technology Conference, Houston, TX; 2010. (13) Moridis, G. J. Numerical Studies of Gas Production From Class 2 and Class 3 Hydrate Accumulations at the Mallik Site, Mackenzie Delta, Canada. SPE Reservoir Eval. Eng. 2004, 7, 175−183. (14) Moridis, G. J.; Reagan, M. T.; Boswell, R.; Collett, T. S.; Zhang, K. Preliminary Evaluation of the Production Potential of Recently Discovered Hydrate Deposits in the Gulf of Mexico. Offshore Technology Conference, Houston, TX; 2010. (15) Moridis, G. J.; Reagan, M. T.; Kim, S. J.; Seol, Y.; Zhang, K. Evaluation of the Gas Production Potential of Marine Hydrate Deposits in the Ulleung Basin of the Korean East Sea. SPE J. 2009, 14, 759−781. (16) Moridis, G. J.; Silpngarmlert, S.; Reagan, M. T.; Collett, T.; Zhang, K. Gas Production from a Cold, Stratigraphically Bounded Hydrate Deposit at the Mount Elbert Site, North Slope, Alaska. Mar. Pet. Geol. 2011, 28, 517−534. (17) Myshakin, E. M.; Gaddipati, M.; Rose, K.; Anderson, B. J. Numerical Simulations of Depressurization-Induced Gas Production from Gas Hydrate Reservoirs at the Walker Ridge 313 Site, Northern Gulf of Mexico. Mar. Pet. Geol. 2012, 34, 169−185. (18) Uddin, M.; Wright, F.; Coombe, D. A. Numerical Study of Gas Evolution and Transport Behaviors in Natural Gas−Hydrate Reservoirs. J. Can. Pet. Technol. 2011, 50, 70−89. (19) Wilson, S. J.; Hunter, R. B.; Collett, T. S.; Hancock, S.; Boswell, R.; Anderson, B. J. Alaska North Slope Regional Gas Hydrate Production Modeling Forecasts. Mar. Pet. Geol. 2011, 28, 460−477. (20) Anderson, B.; Hancock, S.; Wilson, S.; Enger, C.; Collett, T.; Boswell, R.; Hunter, R. Formation Pressure Testing at the Mount Elbert Gas Hydrate Stratigraphic Test Well, Alaska North Slope: Operational Summary, History Matching, and Interpretations. Mar. Pet. Geol. 2011, 28, 478−492. (21) Anderson, B. J.; Wilder, J. W.; Kurihara, M.; White, M. D.; Moridis, G. J.; Wilson, S. J.; Pooladi-Darvish, M.; Masuda, Y.; Collett, T. S.; Hunter, R. B.; Narita, H.; Rose, K.; Boswell, R. Analysis of Modular Dynamic Formation Test Results from the Mount Elbert-01 Stratigraphic Test Well, Milne Point Unit, North Slope Alaska. 6th International Conference on Gas Hydrates, Vancouver, BC, Canada; 2008. (22) Moridis, G.; Collett, T. S.; Dallimore, S. R.; Inoue, T.; Mroz, T. Analysis and Interpretation of the Thermal Test of Gas Hydrate Dissociation in the JAPEX/JNOC/GSC et al. Mallik 5L-38 Gas Hydrate Production Research Well. Scientific Results from the Mallik 2007 Gas Hydrate Production Research Well Program, Malkenzie Delta, Northwest Territories, Canada; Geological Survey of Canada: Ottowa, ON, 2005. (23) Bai, Y.; Li, Q.; Zhao, Y.; Li, X.; Du, Y. The Experimental and Numerical Studies on Gas Production from Hydrate Reservoir by Depressurization. Transp. Porous Media 2009, 443−468. (24) Gupta, A.; Moridis, G. J.; Kneafsey, T. J.; Sloan, E. D. Modeling Pure Methane Hydrate Dissociation Using a Numerical Simulator from a Novel Combination of X-ray Computed Tomography and Macroscopic Data. Energy Fuels 2009, 23, 5958−5965. (25) Konno, Y.; Oyama, H.; Nagao, J.; Masuda, Y.; Kurihara, M. Numerical Analysis of the Dissociation Experiment of Naturally Occurring Gas Hydrate in Sediment Cores Obtained at the Eastern Nankai Trough, Japan. Energy Fuels 2010, 24, 6353−6358. (26) Li, G.; Li, B.; Li, X.-S.; Zhang, Y.; Wang, Y. Experimental and Numerical Studies on Gas Production from Methane Hydrate in Porous Media by Depressurization in Pilot-Scale Hydrate Simulator. Energy Fuels 2012, 26, 6300−6310. (27) Moridis, G.; Seol, Y.; Kneafsey, T. J. Studies of Reaction Kinetics of Methane Hydrate Dissocation in Porous Media. 5th International Conference on Gas Hydrates, Trondheim, Norway; 2005. (28) Ruan, X.; Song, Y.; Liang, H.; Yang, M.; Dou, B. Numerical Simulation of the Gas Production Behavior of Hydrate Dissociation by 5586
dx.doi.org/10.1021/ef500255y | Energy Fuels 2014, 28, 5573−5586
