Effect of Permeability on Hydrate-Bearing Sediment Productivity and Stability in Ulleung Basin, East Sea, South Korea

Methane hydrate has attracted attention as a next-generation resource, and many researchers have conducted various studies to estimate its productivity. Numerical simulation is the optimal method for estimating methane gas productivity. Meanwhile, using a reasonable input parameter is essential for obtaining accurate numerical modeling results. Permeability is a geotechnical property that exhibits the greatest impact on productivity. The permeability of hydrate-bearing sediment varies based on the sediment pore structure and hydrate saturation. In this study, an empirical permeability model was derived from experimental data using soil specimens from the Ulleung Basin, and the model was applied in numerical analysis to evaluate the sediment gas productivity and ground stability. The gas productivity and stability of hydrate-bearing sediments were compared by applying a widely used permeability model and the proposed model to a numerical model. Additionally, a parametric study was performed to examine the effects of initial hydrate saturation on the sediment gas productivity and stability. There were significant differences in the productivity and stability analysis results according to the proposed permeability model. Therefore, it was found that for accurate numerical analysis, a regional permeability model should be applied.


Introduction
Clathrate hydrates are compounds of guest molecules, such as methane, ethane, and propane, with water, and they form at high pressure (typically over 0.6 MPa) and low temperature (typically less than 300 K). Therefore, hydrates are readily formed in the deep ocean, which exhibits high pressure and low temperatures. The guest molecules (also called natural gas) in clathrate hydrates are considered to be premium fuel because they burn cleanly and produce less CO 2 [1]. Also, gas hydrate is considered to be one of the new clean energy resources, which can replace coal, petroleum and conventional natural gas in the 21st century [2]. According to previous research, approximately 230 natural gas hydrate deposits have been investigated globally, with reserves of approximately 1.5 × 10 15 m 3 of natural gas [3]. The commercial production of 15% of the gas hydrate would provide the world with energy for 200 years at the current energy consumption level. The production of natural gas from the hydrates could contribute not only to sustained economic development of individual countries, but also to the global warming mitigation [2]. Various methodologies have been proposed to extract these natural gases from hydrate-bearing sediments, including heat, inhibitor injection, depressurization, CO 2 /CH 4 exchange, and combined methods [4,5]. However, the field-scale production Table 1. Effective permeability models.

Theoretical Model Parameter Remark
Reference Grain coating Kleinberg [20,21] Pore filling Capillary tube Masuda [22] Note: K eff is the effective permeability in the presences of hydrate, K int is the intrinsic permeability of sediment without hydrate, S h is the hydrate saturation.
In this study, a permeability test was conducted using an artificial sample that imitated Ulleung basin sediments. The permeability test modeled the hydrate production process. The permeability data were compared with those of existing models (e.g., Kleinberg and Masuda models) and were used to derive the experimental model for the Ulleung Basin. Comparative numerical analyses using two-dimensional fast Lagrangian analysis of continua (FLAC2D) were performed to estimate the gas and water productivities and stability of the sediments according to the permeability models and initial hydrate satura- tion. The productivity and stability results could be used as the reference data in the step of establishing a strategy for pilot filed-production test in the Ulleung basin, South Korea.

Gas Hydrates in Ulleung Basin
Gas hydrate research in Ulleung Basin was initially conducted by KIGAM in 1996 [8].
In 1998, the potential presence of gas hydrates was first identified in Ulleung Basin, South Korea [23]. In 2007, the first Ulleung Basin Gas Hydrates Drilling Expedition (UBGH1) was performed, and the core along with the downhole log data at the depth range from 55.7 to 231.7 m below sea floor (mbsf) was collected. Gas hydrates were found at three core sites (e.g., UBGH1-04, 09 and 10). In the UBGH1-09 and UBGH1-10 sites, vein-filling gas hydrates were formed in the muddy sediments, while pore-filling gas hydrates were found within the turbidite sand layer in UBGH1-04 [24]. The second Ulleung Basin Gas Hydrates Drilling Expedition (UBGH2) was conducted to better identify the overall distribution of gas hydrate-bearing sediments and to find the potential site for future production tests. During the UBGH2 expedition, 13 sites were explored by logging-while-drilling and coring and characterized for their geological, geophysical, geochemical, geomechanical properties. From the acquired data sets, the geothermal gradient of the Ulleung Basin ranges from 97 to 115 • C/km. The average gas hydrate saturations of UBGH2-2_1 and UBGH2-2_2, which are fracture-dominated reservoirs (i.e., chimney structure), are estimated from the isotropic analysis of resistivity and acoustic velocity logs, which were found to be 46% and 45%, respectively. Meanwhile, the gas hydrate saturation of the UBGH2-10, in which the gas hydrates formed as pore-filling type within the turbidite sand reservoirs, was estimated to range between 35% and 57% [25]. Moreover, expeditions UBGH1 and UBGH2 showed that the gas hydrate in the Ulleung Basin was predominantly CH4, which was generated by microbial activity via CO 2 reduction, with a structure I hydrate [26,27].

UBGH2-6
The UBGH2-6 is located at the most northern site of UBGH2 (N 36 • 01 07.69" E 130 • 15 51.48"), and is considered as a feasible pilot production site ( Figure 1) [8]. The water depth of the UBGH2-6 is 2153 m. The UBGH2-6 contains "high-quality" hydrate-bearing sand reservoirs, and the total cumulative thicknesses of the hydrate-bearing sediments is the greatest of those encountered during the UBGH1 and UBGH2 expeditions. The gas hydrate zones (GHZ) of the UBGH2-6 site occur within the 112-154 m depths. The turbidite sands of UBGH2-6 usually consist of well-distributed coarse silt to fine sands. The gas hydrate saturations within the sand-rich reservoirs at the UBGH2-6 site range between 12% and 79%, with an average of approximately 52% [28].

Soil
Due to the scarcity of in situ core samples, artificial samples that imitated core samples obtained from the silty sand layer in UBGH2-6 were used. The artificial specimen was composed of HAMA #8 sand, in which the iron contents present in the sand were removed

Soil
Due to the scarcity of in situ core samples, artificial samples that imitated core samples obtained from the silty sand layer in UBGH2-6 were used. The artificial specimen was composed of HAMA #8 sand, in which the iron contents present in the sand were removed before testing. The artificial sample closely resembled the core sample in particle size distribution, as shown in Figure 2. The sample exhibited a specific gravity of 2.62 and an average particle size of approximately 130 µm.

Experimental Setup
The experimental setup was designed to model high-pressure and low-temperature conditions under which hydrate can reliably form. Although the in situ gas hydrates are predominantly composed of CH4, CO2 gas was used for the experimental as CH4 is an explosive gas and CO2 and CH4 are both type I structure gas hydrates with similar mechanical properties [29]. In addition, as the purpose of the experiment was to measure the effective permeability, in which the solubility difference has negligible influence, the use

Experimental Setup
The experimental setup was designed to model high-pressure and low-temperature conditions under which hydrate can reliably form. Although the in situ gas hydrates are predominantly composed of CH 4 , CO 2 gas was used for the experimental as CH 4 is an explosive gas and CO 2 and CH 4 are both type I structure gas hydrates with similar mechanical properties [29]. In addition, as the purpose of the experiment was to measure the effective permeability, in which the solubility difference has negligible influence, the use of CO 2 was considered appropriate. The experimental setup used in this study is illustrated in Figure 3. High-pressure conditions were maintained using an oedometer cell, and a low temperature was preserved by using a water bath. To measure the sample permeability, a pressure difference was generated between the two cross-sections of the specimen using a CO 2 and water injection pump and a back-pressure regulator. This pressure difference allowed the pressure to be measured in real time by placing a pressure transducer (accuracy: 0.25% BFSL, Best Fit Straight Line) in front of the pump and using a back-pressure regulator. CO 2 gas was injected from the bottom to the top of the sample. A K-type thermocouple (standard limits of error: 2.2 • C or 0.75% above 0 • C, 2.2 • C or 2.0% below 0 • C) was installed in the specimen to measure the specimen temperature. A CO 2 gas-capturing system was also installed: The discharged CO 2 gas volume from the dissociated gas hydrate was directly measured using the CO 2 gas capturing system. of CO2 was considered appropriate. The experimental setup used in this study is illustrated in Figure 3. High-pressure conditions were maintained using an oedometer cell, and a low temperature was preserved by using a water bath. To measure the sample permeability, a pressure difference was generated between the two cross-sections of the specimen using a CO2 and water injection pump and a back-pressure regulator. This pressure difference allowed the pressure to be measured in real time by placing a pressure transducer (accuracy: 0.25% BFSL, Best Fit Straight Line) in front of the pump and using a backpressure regulator. CO2 gas was injected from the bottom to the top of the sample. A Ktype thermocouple (standard limits of error: 2.2 °C or 0.75% above 0 °C, 2.2 °C or 2.0% below 0 °C) was installed in the specimen to measure the specimen temperature. A CO2 gas-capturing system was also installed: The discharged CO2 gas volume from the dissociated gas hydrate was directly measured using the CO2 gas capturing system.

Experimental Procedure
The experiment consisted of measuring the permeability by sequentially decreasing the hydrate saturation. After dissociating a portion of the gas hydrate, a permeability measurement test was conducted and repeated throughout the dissociation stages.

Experimental Procedure
The experiment consisted of measuring the permeability by sequentially decreasing the hydrate saturation. After dissociating a portion of the gas hydrate, a permeability measurement test was conducted and repeated throughout the dissociation stages.

Hydrate Formation
To simulate the conditions of gas hydrate-bearing sediments in the Ulleung Basin, a pneumatic machine was used to apply 600 kPa of pressure as the effective stress that is normally caused by seabed sediments. In this study, experiments were performed at a temperature of~2 • C and pressure of~3 MPa. At a temperature of approximately 2 • C, the gas hydrate is stable at a pressure of 1.7 MPa or greater; therefore, the gas hydrate formation conditions are satisfied. The Ulleung Basin gas hydrate-bearing sediments exhibit high gas hydrate saturations of~60% at the UBGH2-6 site. It is difficult to achieve over 60% gas hydrate saturation in an experimental environment [30]. Therefore, gas hydrate formation was repeated four times to accomplish the formation of gas hydrates with a high degree of saturation.
After the primary gas hydrate was formed, N 2 gas was injected into the sample to remove the remaining CO 2 gas in the pores. Water was then injected at 2-3 times the volume of the voids to saturate the sample, and then additional N 2 gas was injected at a level of twice the pore volume to remove the water from the pores. Water remained in the particle meniscus, and secondary gas hydrate was formed with the remaining water. N 2 gas was injected to remove CO 2 gas and water because the CO 2 gas solubility is high, and the N 2 gas reactivity is low. This process was repeated three more times to achieve a high gas hydrate saturation.

Permeability Measurement
Constant pressure was maintained in the pump that injected the water; at the surface where the water was drained, the pressure inside the test cell was maintained using a back-pressure regulator. The pressure and temperature of the system were continuously measured during entire process in order to monitor the stability of the gas hydrates. Testing was conducted while observing that no CO 2 was released during the experimental procedure, thereby ensuring the stability of the gas hydrates. The water injection pressure was kept constant, and the pressure on the side where the water was drained was lowered to generate a pressure difference about 50 kPa between the two cross-sections; water then flowed along this pressure difference. The system pressure was maintained, and permeability was measured under steady-state conditions. The amount of water drained after back-pressure regulation and the flow rate through the sample per unit of time were measured. Based on these results, the permeability coefficient was calculated according to the Darcy law [31]. The permeability measurement was repeated for three different porosity values to see the effect of porosity on permeability.

Hydrate Dissociation
The gas hydrate was dissociated by depressurizing the test cell and controlling the back-pressure regulator. The pressure of the system was decreased until the dissociation of the hydrate specimens was observed. After the target volume of the gas was captured, the pressure of the test cell was increased at a constant temperature. Permeability was measured with a 50 kPa head difference at each step of gas hydrate dissociation. After each step, the samples were stabilized for 24 h to ensure that all gases in the pores flowed out and saturated with water. Testing was performed after each stabilization step. Gas hydrate dissociation occurred in 10 steps from the degree of hydrate saturation (S h ) = 0.6 to 0.0. The amount of hydrate dissociated in each step was calculated by collecting the generated CO 2 gas. After gas hydrate dissociation was completed, the initial gas hydrate saturation was estimated using the total volume of the captured CO 2 gas. The experimental process was conducted for 100 h. The water flow rate was low; the CO 2 gas was therefore assumed to be Energies 2021, 14, 1752 6 of 16 fully dissolved in the water. The molar fraction of CO 2 in the water under the experimental conditions was approximately 0.016 [32], which is approximately 12 times higher than the solubility of CO 2 under atmospheric conditions. Therefore, using the captured volume of CO 2 gas may overestimate the specimen gas hydrate saturation. This solubility difference was considered for the final gas hydrate saturation estimation.

Experimental Results and Discussion
The permeability results for various hydrate saturations are displayed in Figure 4a, where K int is the permeability at no hydrate saturation and K eff is the permeability at a certain hydrate saturated condition. During the initial dissociation stage, from 0.6, the permeability did not change significantly, which is presumed to be due to the absence of local dissociation or an interstitial pore network. In contrast, when the gas hydrate saturation reduced to 0.3-0.4, the permeability sharply increased. A small change of porosity does not significantly affect the tendency of permeability according to the hydrate saturation. Based on these results, the permeability of the hydrate-bearing sediments was greatly affected by the hydrate saturation.
The results of the permeability measurement test performed in this study (porosity = 0.48) were compared with results found in the literature (Figure 4b), and the measured permeability data is summarized in Table 2. Comparison between the permeability data was conducted with the Mallik 2L-38 sand sample and the silica sand No. 7 sample [33]. Even though the intrinsic permeability is different with each soil specimen, the permeability reduction trends with increasing hydrate saturation are similar for all samples.   As summarized in Table 1, there are both theoretical and empirical models for permeability and gas hydrate saturation that are widely used for gas productivity numerical modeling. Figure 5 compares the normalized permeability measured by the experiment with that estimated by existing models (i.e., Kleinberg and Masuda models). In contrast to experimental results, the existing models display a tendency for the water permeability to continuously decrease, even up to 70% hydrate saturation. This decrease in permeability as the hydrate saturation increases is a result of the reduction in the cross-sectional area of the pore channel, where water and gas flow, caused by the solid hydrate. The existing theoretical models assumed that the hydrates are formed along the wall (Masuda model) or within the pore center (   As summarized in Table 1, there are both theoretical and empirical models for permeability and gas hydrate saturation that are widely used for gas productivity numerical modeling. Figure 5 compares the normalized permeability measured by the experiment with that estimated by existing models (i.e., Kleinberg and Masuda models). In contrast to experimental results, the existing models display a tendency for the water permeability to continuously decrease, even up to 70% hydrate saturation. This decrease in permeability as the hydrate saturation increases is a result of the reduction in the cross-sectional area of the pore channel, where water and gas flow, caused by the solid hydrate. The existing theoretical models assumed that the hydrates are formed along the wall (Masuda model) or within the pore center (Kleinberg model) in the capillary tube. In addition, both existing models assumed the laminar regime [21]. However, the realistic flow conditions would be closer to turbulence flow conditions. Typically, the permeability in the turbulence flow is lower than that of laminar flow at the same pore space. Moreover, silty sand has a small pore size, and low intrinsic permeability. This leads to a significant decrease in the permeability as hydrates form in the pore, as smaller pore sizes are more largely affected by changes in the pore space. As such, inconsistencies in the experimental and theoretical models may be reduced should these parameters be considered in the existing models. Moreover, the existing models assume that the hydrates are impervious (K eff = 0); however, should a realistic value of permeability be used, it is assumed that a similar converging behavior will be observed with the existing models. In UBGH2-6, the sediments consist of silty sand, which has a relatively lower intrinsic permeability. Furthermore, the impervious permeability of silt is typically considered as 1 mD [34]. However, the existing theoretical models do not take into account the impervious permeability of silty sands.

Permeability Model Derivation
As shown in Figure 5, there was a significant difference between the permeability trends derived from the existing models and those obtained from the experiment using the Ulleung Basin sample. In the case of UBGH2-6, the initial hydrate saturation is approximately 65%. Considering that water permeability is the most critical geotechnical property for gas productivity, the gas productivity can be expected to be overestimated when applying the existing models to a numerical analysis of gas productivity of hydratebearing sediments in the Ulleung Basin. To obtain accurate gas productivity estimations, a permeability model was derived from the experimental results as follows: where

Permeability Model Derivation
As shown in Figure 5, there was a significant difference between the permeability trends derived from the existing models and those obtained from the experiment using the Ulleung Basin sample. In the case of UBGH2-6, the initial hydrate saturation is approximately 65%. Considering that water permeability is the most critical geotechnical property for gas productivity, the gas productivity can be expected to be overestimated when applying the existing models to a numerical analysis of gas productivity of hydratebearing sediments in the Ulleung Basin. To obtain accurate gas productivity estimations, a permeability model was derived from the experimental results as follows: K e f f = K int (0.0045 + 0.9955 where K eff is the effective permeability [mD], K int is the intrinsic permeability [mD] (herein, K int is measured as 209 mD), and S h is the hydrate saturation [-]. where Keff is the effective permeability [mD], Kint is the intrinsic permeability [mD] (herein, Kint is measured as 209 mD), and Sh is the hydrate saturation [-]. Figure 6 shows application of the proposed experimental model on other experimental test data (silica sand No. 7 and Mallik 2L-28). The proposed experimental model is based on the logistic function. The proposed experimental model well matched other experimental data by adjusting the model parameters. This result shows that the proposed experimental model is applicable to other hydrate-bearing sediments if the necessary permeability data are available.

Numerical Simulator
The scale of the gas production mechanism of hydrate-bearing sediments is massive, making it difficult to estimate the gas productivity and stability of the sediments through an experiment. Therefore, it is essential to perform numerical analysis to analyze the gas production mechanism. In this study, a two-dimensional thermo-hydro-mechanical (THM) simulator developed by Kim et al. [35,36] was used to estimate the sediment gas productivity and stability. The simulator was developed based on FLAC2D. FLAC2D is a widely used commercial program in the geotechnical field. FLAC2D cannot perform multi-phase flow coupled with thermal analysis using its built-in constitutive models [37]. However, the program offers a user-defined code module called FISH, which can greatly expand the applicability of FLAC according to user knowledge.

Model Conditions
The analysis target area modeled the UBGH2-6 site, one of the proposed test production sites in the Ulleung Basin in South Korea. The hydrate-bearing sediments are located between the unconsolidated sand layers in the depth range of 140-153 mbsf under the seawater of 2156 m. The detailed layer geometry is illustrated in Figure 7. The permeability models used the experimentally determined model (Equation (1)) and the Kleinberg model [21] for the gas productivity comparison. The Kleinberg model was selected as the control group because a prior core analysis of the Ulleung Basin confirmed that most of the hydrates in the sand layer were distributed in the form of pore filling [38]. The intrinsic permeability of the sand sediment was designated as 209 mD. The initial hydrate saturation was set to 65%. Due to a lack of experiments on stiffness and hydrate saturation for the Ulleung Basin, the Uchida model was used to consider the stiffness change according to the hydrate saturation [39]. The gas production method was modeled using a depressurization method. Depressurization was conducted in a depth range of 140-153 mbsf, and the depressurization rate was 0.5 Mpa/h until the bottom hole pressure was 9 MPa. The effects of sand production were ignored in this study. The production period was set to 14 days. The detailed input properties used in the numerical analysis are summarized in Table 3.
Energies 2021, 14, 1752 9 of 16 tion was set to 65%. Due to a lack of experiments on stiffness and hydrate saturation for the Ulleung Basin, the Uchida model was used to consider the stiffness change according to the hydrate saturation [39]. The gas production method was modeled using a depressurization method. Depressurization was conducted in a depth range of 140-153 mbsf, and the depressurization rate was 0.5 Mpa/h until the bottom hole pressure was 9 MPa. The effects of sand production were ignored in this study. The production period was set to 14 days. The detailed input properties used in the numerical analysis are summarized in Table 3.

Effects of Permeability Models on Water and Gas Productivity
The estimates of water pumped out and gas production according to the permeability models (e.g., experimental and Kleinberg models) are illustrated in Figure 8. In the Kleinberg model, the early stages of production (up to 1 day) exhibited a rapid increase in the water production rate, followed by a rapid decline. In contrast, for the experimental model, the water production rate exhibited a relatively smooth rise. The experimental model (S h = 30% and 65%) exhibited a low initial permeability of the hydrate-bearing sediments compared to that of the Kleinberg model (S h = 65%) ( Figure 5). Initially the water and gas production rate likely increased rapidly due to high permeability and afterward rapidly decreased due to ground subsidence. In other words, the initial stage is the state in which the pore pressure continuously decreases towards the target bottom hole pressure. The pressure gradient is relatively high and causes an increase the specific water discharge. An increase in effective stress, caused by the decrease in pore water pressure, induced ground subsidence, and it caused the decrease in permeability at the end of the depressurization. models (e.g., experimental and Kleinberg models) are illustrated in Figure 8. In the Kleinberg model, the early stages of production (up to 1 day) exhibited a rapid increase in the water production rate, followed by a rapid decline. In contrast, for the experimental model, the water production rate exhibited a relatively smooth rise. The experimental model (Sh = 30% and 65%) exhibited a low initial permeability of the hydrate-bearing sediments compared to that of the Kleinberg model (Sh = 65%) ( Figure 5). Initially the water and gas production rate likely increased rapidly due to high permeability and afterward rapidly decreased due to ground subsidence. In other words, the initial stage is the state in which the pore pressure continuously decreases towards the target bottom hole pressure. The pressure gradient is relatively high and causes an increase the specific water discharge. An increase in effective stress, caused by the decrease in pore water pressure, induced ground subsidence, and it caused the decrease in permeability at the end of the depressurization.
The maximum water and gas production rates in the Kleinberg model were 12.32 ton/day and 2.76 × 10 3 m 3 /day, respectively. In the experimental model, the maximum water and gas production rates were 6.25 ton/day and 1.06 × 10 3 m 3 /day, respectively. In the experimental model, the predicted production rates for both water and gas were estimated to be approximately 50% and 60% less than those of the Kleinberg model, respectively (Figure 8a,c). The cumulative water and gas production were 115.43 tons and 2.06 × 10 4 m 3 in the Kleinberg model and 72.67 tons and 9.68 × 10 3 m 3 in the experimental model, respectively. Consequently, the experimental model predicted the total water and gas production as approximately 37% and 53% less, respectively, than that of the Kleinberg model (Figure 8b,d). In regards to the numerical results, it shows that the existing model overestimates the productivity in the case of Ulleung Basin.  (c) (d) Figure 8. Water and gas production after 14 days for studied permeability models: (a) water production rate; (b) cumulative water production; (c) gas production rate; (d) cumulative gas production.

Effects of Initial Hydrate Saturation on Productivity
The experimental results confirmed that the permeability of hydrate-bearing sediments depends on the hydrate saturation. Thus, a parametric analysis was performed to investigate the effects of the initial hydrate saturation on productivity. The initial hydrate saturations were set to 30%, 55% and 65% to investigate the effects of the initial hydrate saturation on the sediment gas productivity and stability. The parametric analysis productivity predictions are displayed in Figure 9. The results for an initial hydrate saturation of 65% are the same as those of the experimental model in the previous section. The maximum water and gas production rates for an initial hydrate saturation of 30% were 75.32 ton/day and 2.37 × 10 4 m 3 /day, respectively. The estimated water and gas production rates were approximately 12 times and 22 times more for an initial hydrate saturation of 30% than for that of 65%, respectively (Figure 9a,c). The cumulative water and gas production   Figure 8. Water and gas production after 14 days for studied permeability models: (a) water production rate; (b) cumulative water production; (c) gas production rate; (d) cumulative gas production.
The maximum water and gas production rates in the Kleinberg model were 12.32 ton/day and 2.76 × 10 3 m 3 /day, respectively. In the experimental model, the maximum water and gas production rates were 6.25 ton/day and 1.06 × 10 3 m 3 /day, respectively. In the experimental model, the predicted production rates for both water and gas were estimated to be approximately 50% and 60% less than those of the Kleinberg model, respectively (Figure 8a,c). The cumulative water and gas production were 115.43 tons and 2.06 × 10 4 m 3 in the Kleinberg model and 72.67 tons and 9.68 × 10 3 m 3 in the experimental model, respectively. Consequently, the experimental model predicted the total water and gas production as approximately 37% and 53% less, respectively, than that of the Kleinberg model (Figure 8b,d). In regards to the numerical results, it shows that the existing model overestimates the productivity in the case of Ulleung Basin.

Effects of Initial Hydrate Saturation on Productivity
The experimental results confirmed that the permeability of hydrate-bearing sediments depends on the hydrate saturation. Thus, a parametric analysis was performed to investigate the effects of the initial hydrate saturation on productivity. The initial hydrate saturations were set to 30%, 55% and 65% to investigate the effects of the initial hydrate saturation on the sediment gas productivity and stability. The parametric analysis productivity predictions are displayed in Figure 9. The results for an initial hydrate saturation of 65% are the same as those of the experimental model in the previous section. The maximum water and gas production rates for an initial hydrate saturation of 30% were 75.32 ton/day and 2.37 × 10 4 m 3 /day, respectively. The estimated water and gas production rates were approximately 12 times and 22 times more for an initial hydrate saturation of 30% than for that of 65%, respectively (Figure 9a,c). The cumulative water and gas production for an initial hydrate saturation of 30% were 578.34 tons and 1.72 × 10 5 m 3 , respectively. When the initial hydrate saturation was decreased by about half, the total water and gas production predictions increased by 7.9 times and 17.7 times, respectively (Figure 9b,d). This significant difference is inferred to be due to the difference in the initial permeability and the tendencies of permeability change according to the hydrate saturation. The initial permeability of experimental model (S h = 30%) exhibited a similar value with the Kleinberg model (S h = 0.65%), but the experimental model exhibited a more rapid increase in permeability as the hydrate dissociation progressed than that of the Kleinberg model ( Figure 5). This leads to the significant difference in productivity with production period. (c) (d) Figure 9. Water and gas production after 14 days for different initial gas hydrate saturations where the experimentally determined permeability model was used: (a) water production rate; (b) cumulative water production; (c) gas production rate; (d) cumulative gas production.
In addition, the predicted gas production estimate was indirectly compared with the Mallik data [10]. The initial hydrate saturation in Mallik was high, ranging from 55% to 85%. For more accurate comparison, the parametric case (initial hydrate saturation = 55%) was added while applying the experimental model and the Kleinberg model. The intrinsic permeability of sediment in Mallik was 10 mD to 1000 mD. The average permeability, 500   Figure 9. Water and gas production after 14 days for different initial gas hydrate saturations where the experimentally determined permeability model was used: (a) water production rate; (b) cumulative water production; (c) gas production rate; (d) cumulative gas production. In addition, the predicted gas production estimate was indirectly compared with the Mallik data [10]. The initial hydrate saturation in Mallik was high, ranging from 55% to 85%. For more accurate comparison, the parametric case (initial hydrate saturation = 55%) was added while applying the experimental model and the Kleinberg model. The intrinsic permeability of sediment in Mallik was 10 mD to 1000 mD. The average permeability, 500 mD, was applied to the comparative analysis. The gas production rates in Mallik 2007 and 2008 were 1000-2000 m 3 /day for 1.5 days and 2000-3000 m 3 /day for 6 days, respectively. The gas production value in Mallik were intermediate to the Kleinberg model (S h = 0.55, K int = 500 mD) and the experimental model (S h = 0.55, K int = 500 mD), as shown in Figure 10. The results of both methods differ from the actual data, but this is due to the lack of available detailed model information from the Mallik site. If accurate information on the site to be analyzed is not available, the productivity of gas and water can be indirectly estimated through the numerical analysis by applying the Kleinberg model and the proposed experimental model.

Stability of Hydrate-Bearing Sediments
Stability analysis of the hydrate-bearing sediments was conducted according to each permeability model and initial hydrate saturation. During hydrate production, displacement data were collected from three locations: on the sea floor and on the top and bottom of the production well ( Figure 11). In all cases, ground subsidence occurred at the sea floor and the production well head (Figure 11a,b). The case of an initial hydrate saturation of 30% and the experimental model exhibited the largest vertical displacements. In contrast, the case with the lowest ground subsidence was that with the experimental model and an initial hydrate saturation of 65%. The stability analysis results ( Figure 11) also relate to productivity. The decrease in the sediment stiffness is due to the dissociation of the solid hydrate in the pore, which causes the sediment subsidence. Thus, the higher the water and gas productivity, the larger the ground displacement, and the faster it progresses.
In particular, ground heave occurred at the bottom of the production well in all cases (Figure 11c). Similar to the subsidence results, the largest ground heave occurred in the case of the experimental model and an initial hydrate saturation of 30%. In this case, the maximum ground uplift was approximately 35.21 cm. In contrast, for the developed model and an initial hydrate saturation of 65%, the smallest ground heave was 3.28 cm. Based on these results, the subsidence on the sea floor and the top of the production well may differ by up to 9.5 times and 10.8 times, respectively, according to the applied permeability model and initial hydrate saturation. The ground heave at the bottom of the production well is induced by the pressure difference between the discharge pressure and pore pressure. Ground heave and subsidence during gas production can cause tensile or compressive stress on the production well. Thus, it is essential that additional analysis on the stability of production wells during gas production be conducted. Gas production rate [m 3 /day]

Stability of Hydrate-Bearing Sediments
Stability analysis of the hydrate-bearing sediments was conducted according to each permeability model and initial hydrate saturation. During hydrate production, displacement data were collected from three locations: on the sea floor and on the top and bottom of the production well ( Figure 11). In all cases, ground subsidence occurred at the sea floor and the production well head (Figure 11a,b). The case of an initial hydrate saturation of 30% and the experimental model exhibited the largest vertical displacements. In contrast, the case with the lowest ground subsidence was that with the experimental model and an initial hydrate saturation of 65%. The stability analysis results ( Figure 11) also relate to productivity. The decrease in the sediment stiffness is due to the dissociation of the solid hydrate in the pore, which causes the sediment subsidence. Thus, the higher the water and gas productivity, the larger the ground displacement, and the faster it progresses.
In particular, ground heave occurred at the bottom of the production well in all cases (Figure 11c). Similar to the subsidence results, the largest ground heave occurred in the case of the experimental model and an initial hydrate saturation of 30%. In this case, the maximum ground uplift was approximately 35.21 cm. In contrast, for the developed model and an initial hydrate saturation of 65%, the smallest ground heave was 3.28 cm. Based on these results, the subsidence on the sea floor and the top of the production well may differ by up to 9.5 times and 10.8 times, respectively, according to the applied permeability model and initial hydrate saturation. The ground heave at the bottom of the production well is induced by the pressure difference between the discharge pressure and pore pressure.
Ground heave and subsidence during gas production can cause tensile or compressive stress on the production well. Thus, it is essential that additional analysis on the stability of production wells during gas production be conducted.

Limitations and Further Work
The experimental model was derived from the lab-scale permeability test. However, the effective stress cannot be increased indefinitely in laboratory tests. The permeability test was conducted at an effective stress of 600 kPa; however, the effective stress in the Ulleung Basin is 14 MPa after depressurization. Therefore, a correlation between the confining stress level and hydraulic properties of the HBS is required for accurate numerical analysis.
The Young's modulus of hydrate-bearing sediments varies based on the hydrate saturation. During hydrate dissociation, the hydrate-bearing sediment stiffness decreases. Therefore, numerical analyses must consider the effects of hydrate saturation on the Young's modulus of sediments. Similar to the permeability model, the stiffness model may also depend on the sediment condition. Therefore, deriving a stiffness model for the Ulleung Basin and applying it to a numerical model is necessary.
Sand production effects were not considered in this study. In the first field production test of the Mallik project in 2007, the production period was only 1.5 days due to the occurrence of sand production [10]. Sand production can lead to production well clogging, and excessive stress concentrated around the production well can cause the well to collapse [43]. Therefore, in future studies, productivity and stability analyses that consider the effects of sand production are required.   Figure 11. Subsidence and heave of sediment near the production well: (a) subsidence on sea floor; (b) subsidence at production well head; (c) ground heave at production well bottom.

Limitations and Further Work
The experimental model was derived from the lab-scale permeability test. However, the effective stress cannot be increased indefinitely in laboratory tests. The permeability test was conducted at an effective stress of 600 kPa; however, the effective stress in the Ulleung Basin is 14 MPa after depressurization. Therefore, a correlation between the confining stress level and hydraulic properties of the HBS is required for accurate numerical analysis.
The Young's modulus of hydrate-bearing sediments varies based on the hydrate saturation. During hydrate dissociation, the hydrate-bearing sediment stiffness decreases. Therefore, numerical analyses must consider the effects of hydrate saturation on the Young's modulus of sediments. Similar to the permeability model, the stiffness model may also depend on the sediment condition. Therefore, deriving a stiffness model for the Ulleung Basin and applying it to a numerical model is necessary.
Sand production effects were not considered in this study. In the first field production test of the Mallik project in 2007, the production period was only 1.5 days due to the occurrence of sand production [10]. Sand production can lead to production well clogging, and excessive stress concentrated around the production well can cause the well to collapse [43]. Therefore, in future studies, productivity and stability analyses that consider the effects of sand production are required.

Conclusions
This study examined the effects of permeability and hydrate saturation on the gas productivity and stability of hydrate-bearing sediments. The permeability test simulated the hydrate production mechanism, and permeability and hydrate saturation were measured. Using the measured data, an experimental model was derived and applied to the numerical analysis to estimate the sediment water and gas productivities and stability. The key findings of this study are as follows:

•
The measured permeability using the artificial specimen of the Ulleung Basin has significant differences compared to that of theoretical models. The measured permeability exhibited differences of up to 75 times in the Masuda model (N = 10) and approximately 15 times in the Kleinberg model (grain coating), at 60% hydrate saturation. • The experimental model showed the convergence in permeability as the hydrate saturation increase, reflecting the impervious permeability of silty sand, whereas the existing theoretical models did not reflect the hydraulic characteristics of the sediments, and showed a continuous decrease in permeability as the hydrate saturation increases.

•
In the numerical estimation that applied the experimental model, the maximum water and gas production were approximately 50% and 60% smaller than those of the Kleinberg model. The cumulative water and gas production for the experimental model predicted the total water and gas production as approximately 37% and 53% less, respectively, than those of the Kleinberg model.

•
The initial hydrate saturation significantly affects the productivity. This is due to the initial permeability and the tendencies of permeability change due to the hydrate saturation. This means that in long-term production, the productivity might increase as the hydrate dissociates from the high initial hydrate saturation, which in certain production periods can show inverse productivity tendencies in regards to the initial hydrate saturation.

•
The hydrate-bearing sediment stability was similar to the production estimates. The ground displacement tended to be proportional to the permeability and inversely proportional to the initial hydrate saturation.

•
Ground heave occurs at the bottom of the production well, the interaction between the ground heave, and subsidence might induce tensile or compressive stresses on the production well.
The existing model over-estimated the water and gas productivities and stability of the hydrate-bearing sediments. This research improves the reliability of the gas productivity and stability predictions for the Ulleung Basin. However, this study also has limitations, such as the effects of sand production not being considered. This means that productivity and stability may be slightly overestimated. Further research is required to develop a numerical model that considers sand production to increase productivity and stability evaluation reliability.

Data Availability Statement:
The data presented in this study are available on request from the corresponding author.