Wave Propagation Characteristics in Gas Hydrate-Bearing Sediments and Estimation of Hydrate Saturation

: Natural gas hydrate is a new clean energy source in the 21st century, which has become a research point of the exploration and development technology. Acoustic well logs are one of the most important assets in gas hydrate studies. In this paper, an improved Carcione–Leclaire model is proposed by introducing the expressions of frame bulk modulus, shear modulus and friction coefﬁcient between solid phases. On this basis, the sensitivities of the velocities and attenuations of the ﬁrst kind of compressional (P1) and shear (S1) waves to relevant physical parameters are explored. In particular, we perform numerical modeling to investigate the effects of frequency, gas hydrate saturation and clay on the phase velocities and attenuations of the above ﬁve waves. The analyses demonstrate that, the velocities and attenuations of P1 and S1 are more sensitive to gas hydrate saturation than other parameters. The larger the gas hydrate saturation, the more reliable P1 velocity. Besides, the attenuations of P1 and S1 are more sensitive than velocity to gas hydrate saturation. Further, P1 and S1 are almost nondispersive while their phase velocities increase with the increase of gas hydrate saturation. The second compressional (P2) and shear (S2) waves and the third kind of compressional wave (P3) are dispersive in the seismic band, and the attenuations of them are signiﬁcant. Moreover, in the case of clay in the solid grain frame, gas hydrate-bearing sediments exhibit lower P1 and S1 velocities. Clay decreases the attenuation of P1, and the attenuations of S1, P2, S2 and P3 exhibit little effect on clay content. We compared the velocity of P1 predicted by the model with the well log data from the Ocean Drilling Program (ODP) Leg 164 Site 995B to verify the applicability of the model. The results of the model agree well with the well log data. Finally, we estimate the hydrate layer at ODP Leg 204 Site 1247B is about 100–130 m below the seaﬂoor, the saturation is between 0–27%, and the average saturation is 7.2%.


Introduction
Natural gas hydrates are recognized as a kind of potential energy source [1][2][3][4]. However, gas hydrate resources are highly dispersed, but some individual gas hydrate accumulations may hold important and concentrated resources that may be economically exploited in the future [5,6]. Gas hydrates are extremely unstable and easy to decompose during the exploitation process. As a result, the marine sediment structure can be destroyed, leading to geological disasters such as blowouts and collapses, which pollute the marine ecological environment [7][8][9][10]. In addition, methane produced by hydrate decomposition can cause a serious greenhouse effect [11][12][13][14]. Consequently, the exploration and development technology of natural gas hydrate is a research focus and hotspot [15][16][17][18][19]. During the exploitation of this unconventional and clean energy, it is significant to identify the position of the gas hydrate layer and to estimate the hydrate saturation in the near future, which requires the development of appropriate detection technology [20]. Drilling, coring, well logging are three commonly used tools in exploration for natural gas hydrates. However, drilling and coring are things of great difficulty due to the mentioned severe conditions. Hence, well logging is considered to be the priority method for effective exploration and quantitative description of gas hydrate-bearing sediments [21].
As one of the most effective gas hydrate detection methods, the acoustic well log has been studied by many researchers [22][23][24][25][26][27][28][29][30]. It can provide estimates of gas hydrate saturation in a continuous depth profile by using predictive models that relate the hydrate saturation to various parameters such as P-and S-wave velocities or attenuations.
Many models have been proposed for calculating the wave velocities in gas hydratebearing sediments such as Lee's weight equation [31], K-T equation [32], cementation model [33], effective medium theory [34], Biot-Gassmann model [35] and Lee's threephase Biot equation [36]. Much work has been done to study the characteristics of gas hydrate-bearing sediments by using various models mentioned above. Priest et al. [37] used acoustic waves to test the distribution and content of hydrate in gas hydrate-bearing sediments, and analyzed the effect of hydrate saturation on velocity and attenuation amplitude. Zhang [38] and Zhang [39], etc., equalized the three-phase porous medium as a Biot two-phase medium, and studied the relationship between three-wave velocities and attenuations with frequency.
The above models to estimate hydrate saturation have their strengths and limitations. Though these models can estimate hydrate saturation to a certain extent, all of them simplify hydrate reservoirs that don't agree with the real situation. To precisely predict the distribution of the hydrate and estimate the hydrate saturation, it's necessary to establish a predictive petrophysical model of describing the acoustic characteristics of the gas hydrate. In fact, gas hydrate sediment is a typical three-phase porous medium, which is composed of a solid grain frame with a fluid (usually is water) and a solid coexist in pores. Leclaire [40] proposed the percolation theory based on Biot's theory [41,42] and analyzed the wave propagation in frozen porous media (frozen soil or permafrost). However, Leclaire assumed that there was no direct contact between solid grains and ice. Carcione [43] proposed the Carcione-Leclaire model, considering the contact between two solids, and used the K-T theory to calculate the frame modulus of the ice. Both Carcione and Leclaire pointed out when waves propagate in porous media with two solid and one liquid, three compressional waves (one fast P wave and two slow P waves) and two shear waves (one fast S wave and one slow S wave) can be excited. Carcione-Leclaire model derives the phase velocities and attenuations of five waves in three-phase porous media, which explains the wave propagation better than other models. In practice, slow waves are difficult to be observed, but they slack off a part of the energy, which should not be ignored. Carcione neglects friction between two solids, which limits its practical application. Geurin et al. [44] defined the friction coefficient between two solids based on the drag coefficient calculated by Berryman and Wang [45]. Besides, Lee et al. [46] proposed a more straightforward and more effective method to calculate skeleton modulus by using two free parameters. The completeness of Carcione-Leclaire model was perfected by adding the expressions of these parameters into the model.
Herein, we propose an improved Carcione-Leclaire model aiming to analyze the wave characteristics in gas hydrate-bearing sediments and provide saturation estimation. We first discuss the sensitivities of velocities and attenuations of P1 and S1 to various parameters especially gas hydrate saturation. Further, the variation of phase velocities and attenuations of waves with frequency and gas hydrate saturation is emphatically studied. Then the effects of clay on wave propagation characteristics are analyzed in detail. Moreover, we test and verify the model using well log data from ODP Leg 164 Site 995B. We then estimate the hydrate saturation of ODP Leg 204 Site 1247B using the model. Finally, we compare different P-wave velocities of three models (effective medium model, three-phase Biot equation and the improved Carcione-Leclaire model), which is displayed in Appendices A and B.

Improved Carcione-Leclaire Model
Based on Biot's theory [41], the motion equation of three-phase porous media can be derived by using the Lagrange equation, which can be expressed in a matrix form: where , A, R and µ are the matrices of mass densities, friction coefficients, rigidity and shear moduli, respectively, which are given by: where u = (u s i , u w i , u h i ) T denotes the displacement matrix. The superscripts s, w and h refer to solid grain, pore water and gas hydrate, respectively. The subscripts i,j (i, j = 1, 2, 3, i = j) in matrices refer to coupling or interaction between phases i and j. ρ ij are the generalized mass densities. b 12 , b 23 and b 13 represent friction coefficients, 1, 2 and 3 refer to solid grain, pore water and gas hydrate, respectively: where η w is the water viscosity, φ w denotes the volume fraction of the water, κ s and κ h are the effective permeabilities of the solid grain frame and hydrate frame: where φ is the porosity, φ h is the volume fraction of the gas hydrate, φ = φ h + φ w . According to Berryman and Wang [45], the coefficient between solid grains and gas hydrate b 13 is defined as [44]: where the value of the coefficient b 0 13 = 2.2 × 10 8 is determined, which was suggested by Guerin [44].
The matrices of rigidity and shear moduli can be expressed as: where K i and µ i (i = 1, 2, 3) are the bulk modulus and shear modulus of each phase, which are given by: where K sm , K hm , µ sm and µ hm are the bulk modulus and shear modulus of the solid grain frame and hydrate frame, respectively, K av and µ av are average moduli defined by: where c 1 , g 1 , c 3 and g 3 are the consolidation coefficients of solid grains and gas hydrate defined according to the Biot theory: In order to complete the three-phase model, it's necessary to express the elastic moduli of the solid grain frame and gas hydrate frame. Lee et al. [22] express them as: where φ ah and φ as denote the apparent porosities, . ε and α refer to the free parameters. ε accounts for the support of hydrate to the frame. If ε is taken as 0, it implies the hydrate is completely in the form of solid skeleton, and if 1 is taken, it indicates the hydrate is completely suspended in the pores. ε has little influence on the calculation results and we take ε = 0.12, which was suggested by previous studies [46,47]. The consolidation coefficient α indicates the cementation of the whole skeleton, which is a sediment-dependent parameter varying based on the properties of the gas hydrate reservoir, and has a certain effect on the theoretical calculation velocity [48]. α can be calculated using the formula presented by Zhang [49]: where µ ma is the frame shear modulus, ρ is the density of the sediment, V s is the S-wave velocity, S h is the gas hydrate saturation.

Dispersion Equations
Applying the divergence and rotational operators to the motion equation (Equation (1)), the governing equations of compressional and shear waves can be obtained: The displacement potential functions of compressional and shear waves are: By substituting them into Equation (12), the governing equations of compressional and shear waves are: Therefore, we can obtain the characteristic equation of the compressional waves: where Z 1 , Z 2 , Z 3 and Z 4 are given by: where: Similarly, the dispersion equation of shear waves can be obtained: where: The velocities can be expressed as: The attenuation factors can be written as: We use the inverse quality factor Q −1 to evaluate the attenuations of the waves, which can be written as: where j = P1, P2, P3, S1 and S2.

Model Parameters
According to the dispersion equations, there are three compressional and two shear waves propagating in gas hydrate-bearing sediments which can be labeled as P1, P2, P3, S1 and S2, respectively [50]. In this section, we discuss the influence of different parameters on P1 and S1 velocities and attenuations at different frequencies by calculating the sensitivity. The physical parameters of the gas hydrate-bearing sediment constituents are presented in Table 1. According to the description of the consolidation coefficient α in Section 2.1, α is related to the S-wave velocity, which is unknown in the process of the simulation. Hence, we consider a fixed α when calculating the sensitivity and studying the wave propagation characteristics. The phase velocity of P1 as a function of gas hydrate saturation with different α is shown in Figure 1. Two frequencies 100 Hz and 10 kHz are considered. It indicates that the consolidation coefficient α has almost no effect on the trend of P1-velocity variation. Moreover, except for the case α is 20, the P1 velocities obtained under the other three consolidation coefficients are in consistent. Thus, the consolidation coefficient α is set to be 30 in the following numerical simulations [36,37].

Sensitivity Analysis
The reservoir parameters have remarkable effects on the characteristics of the sediment. Sensitivity can be used to describe the influence of parameters on the phase velocity or attenuation at different frequencies. In this section, we discuss the sensitivities of phase velocities and attenuations of P1 and S1. The sensitivity can be expressed as [51]: where v represents the velocity, S denotes the sensitivity, p refers to the medium parameter. The sensitivity is greater or less than zero means that the velocity and attenuation increase or decrease as the variable increases. Figure 2 shows the sensitivities of phase velocities and attenuations of P1 and S1 to hydrate saturation, porosity and three friction coefficients. Porosities of 40-80% fundamentally cover the scope observed in hydrate-bearing sediments [3]. Thus, consider the case that the basic values of the hydrate saturation and the porosity are both chosen to be 0.5. For P1 velocity, Figure 2a indicates that the sensitivity of P1 velocity to porosity and hydrate saturation is higher than that of the three friction coefficients. In particular, the increase of hydrate saturation increases the P1 velocity whereas the increase of porosity decreases the P1 velocity. Moreover, the friction coefficients have almost no effect on the P1 velocity. For P1 attenuation, the attenuation of P1 is the most sensitive to hydrate saturation, whereas less sensitive to friction coefficients especially the friction coefficient between solid grains and hydrate ( Figure 2b). As can be observed in Figure 2c,d, the sensitivity of the velocity and attenuation of S1 to various parameters has the same trend as that of P1. It is worth noting that the sensitivities of phase velocities of P1 and S1 are independent of frequency to all parameters we considered in Figure 2a,c. Obviously, the attenuation of S1 to gas hydrate saturation is less sensitive than that of P1 (Figure 2b   ies 2021, 14, x FOR PEER REVIEW 9 of (c) (d) Figure 2. The sensitivities of velocities and attenuations of P1 and S1 to various parameters. Different parameters are involved: gas hydrate saturation, porosity and three friction coefficients between different phases. 13 b is the friction coefficient between solid grains and gas hydrate; 12 b is the friction coefficient between solid grains and pore water; 23 b is the friction coefficient between pore water and gas hydrate: (a) P1 velocity; (b) P1 attenuation; (c) S1 velocity; (d) S1 attenuation.
These findings suggest that gas hydrate saturation is the main factor affecting the P and S1 velocities in gas hydrate-bearing sediments. Thus, we further investigated the se sitivities of velocities and attenuations of P1 and S1 to gas hydrate saturation as shown Figure 3 with the porosity φ = 0.5. From Figures 3a,c, it can be seen that with the increa of gas hydrate saturation, the sensitivities of velocities of P1 and S1 to gas hydrate satur tion increase, which are independent of frequency. As a consequence, the inversion b velocity is more reliable with higher gas hydrate saturation. As indicated in Figures 3b, the sensitivities of attenuations of P1 and S1 exhibit steep slopes. As the basic value of g hydrate saturation increases, the sensitivities of attenuations of P1 and S1 decrease in th seismic band, but increase in the ultrasonic band at a fixed frequency. Besides, at a fixe saturation, the sensitivities of attenuations of P1 and S1 initially increase rapidly then d crease as frequency increases. Thus, the inversion by attenuation is more reliable in th seismic band especially for lower saturation, while more reliable in the acoustic loggin band for higher saturation.
Sensitivity of P1 velocity Figure 2. The sensitivities of velocities and attenuations of P1 and S1 to various parameters. Different parameters are involved: gas hydrate saturation, porosity and three friction coefficients between different phases. b 13 is the friction coefficient between solid grains and gas hydrate; b 12 is the friction coefficient between solid grains and pore water; b 23 is the friction coefficient between pore water and gas hydrate: (a) P1 velocity; (b) P1 attenuation; (c) S1 velocity; (d) S1 attenuation.
These findings suggest that gas hydrate saturation is the main factor affecting the P1 and S1 velocities in gas hydrate-bearing sediments. Thus, we further investigated the sensitivities of velocities and attenuations of P1 and S1 to gas hydrate saturation as shown in Figure 3 with the porosity φ = 0.5. From Figure 3a,c, it can be seen that with the increase of gas hydrate saturation, the sensitivities of velocities of P1 and S1 to gas hydrate saturation increase, which are independent of frequency. As a consequence, the inversion by velocity is more reliable with higher gas hydrate saturation. As indicated in Figure 3b,d, the sensitivities of attenuations of P1 and S1 exhibit steep slopes. As the basic value of gas hydrate saturation increases, the sensitivities of attenuations of P1 and S1 decrease in the seismic band, but increase in the ultrasonic band at a fixed frequency. Besides, at a fixed saturation, the sensitivities of attenuations of P1 and S1 initially increase rapidly then decrease as frequency increases. Thus, the inversion by attenuation is more reliable in the seismic band especially for lower saturation, while more reliable in the acoustic logging band for higher saturation. (c) (d) Figure 2. The sensitivities of velocities and attenuations of P1 and S1 to various parameters. Different parameters are involved: gas hydrate saturation, porosity and three friction coefficients between different phases. 13 b is the friction coefficient between solid grains and gas hydrate; 12 b is the friction coefficient between solid grains and pore water; 23 b is the friction coefficient between pore water and gas hydrate: (a) P1 velocity; (b) P1 attenuation; (c) S1 velocity; (d) S1 attenuation.
These findings suggest that gas hydrate saturation is the main factor affecting the P1 and S1 velocities in gas hydrate-bearing sediments. Thus, we further investigated the sensitivities of velocities and attenuations of P1 and S1 to gas hydrate saturation as shown in Figure 3 with the porosity φ = 0.5. From Figures 3a,c, it can be seen that with the increase of gas hydrate saturation, the sensitivities of velocities of P1 and S1 to gas hydrate saturation increase, which are independent of frequency. As a consequence, the inversion by velocity is more reliable with higher gas hydrate saturation. As indicated in Figures 3b,d, the sensitivities of attenuations of P1 and S1 exhibit steep slopes. As the basic value of gas hydrate saturation increases, the sensitivities of attenuations of P1 and S1 decrease in the seismic band, but increase in the ultrasonic band at a fixed frequency. Besides, at a fixed saturation, the sensitivities of attenuations of P1 and S1 initially increase rapidly then decrease as frequency increases. Thus, the inversion by attenuation is more reliable in the seismic band especially for lower saturation, while more reliable in the acoustic logging band for higher saturation.  (c) S1 velocity; (d) S1 attenuation.

Numerical Simulation Analysis
Based on the previous discussion, it is obvious that the phase velocities and attenuations of the P1 and S1 are the most sensitive to gas hydrate saturation. In this section, we discuss the phase velocities and attenuations of the five-wave modes versus frequency (especially in the seismic and acoustic logging band) and gas hydrate saturation in detail.

Influence of Frequency and Gas Hydrate Saturation
The calculated velocities of the five modes in gas hydrate-bearing sediments as functions of frequency and gas hydrate saturation are displayed in Figure 4. The sediment porosity in the model is chosen to be 50%. Other parameters are the same as those in the previous section. As indicated in Figures 4a,b, P1 and S1 propagate pretty fast, equivalent to the P and S waves in a single-phase elastic medium. Moreover, there is almost no dispersion of P1 and S1 in the seismic exploration (10-100 Hz) and acoustic logging (1-20 kHz) frequency range, the velocities of these two waves increase almost linearly with the increase of gas hydrate saturation. From Figures 4c,d,e, it can be seen that P2, P3, and S2 propagate slowly especially in the seismic band. This is because in the low-frequency range, the viscosity of pore fluid plays a major role in velocity, while the inertial force has more effect in the high-frequency range. Moreover, the three slow waves are dispersive in the low-frequency range including the seismic band, while there is almost no dispersion in the high-frequency range. At high frequency, the phase velocity of P2 initially decreases and then increases (Figure 4c) as the increase of hydrate saturation, while the phase velocity of S2 increases (Figure 4d), and the phase velocity of P3 increases at first and then decreases (Figure 4e). It is noteworthy that with the increase of hydrate saturation, the transition zones between low frequency and high frequency of P2, S2 and P3 move towards high frequency.

Numerical Simulation Analysis
Based on the previous discussion, it is obvious that the phase velocities and attenuations of the P1 and S1 are the most sensitive to gas hydrate saturation. In this section, we discuss the phase velocities and attenuations of the five-wave modes versus frequency (especially in the seismic and acoustic logging band) and gas hydrate saturation in detail.

Influence of Frequency and Gas Hydrate Saturation
The calculated velocities of the five modes in gas hydrate-bearing sediments as functions of frequency and gas hydrate saturation are displayed in Figure 4. The sediment porosity in the model is chosen to be 50%. Other parameters are the same as those in the previous section. As indicated in Figure 4a,b, P1 and S1 propagate pretty fast, equivalent to the P and S waves in a single-phase elastic medium. Moreover, there is almost no dispersion of P1 and S1 in the seismic exploration (10-100 Hz) and acoustic logging (1-20 kHz) frequency range, the velocities of these two waves increase almost linearly with the increase of gas hydrate saturation. From Figure 4c-e, it can be seen that P2, P3, and S2 propagate slowly especially in the seismic band. This is because in the low-frequency range, the viscosity of pore fluid plays a major role in velocity, while the inertial force has more effect in the high-frequency range. Moreover, the three slow waves are dispersive in the low-frequency range including the seismic band, while there is almost no dispersion in the high-frequency range. At high frequency, the phase velocity of P2 initially decreases and then increases (Figure 4c) as the increase of hydrate saturation, while the phase velocity of S2 increases (Figure 4d), and the phase velocity of P3 increases at first and then decreases (Figure 4e). It is noteworthy that with the increase of hydrate saturation, the transition zones between low frequency and high frequency of P2, S2 and P3 move towards high frequency.
As demonstrated in the previous sections, the attenuation characteristics of acoustic waves can be used to estimate the hydrate saturation. Figure 5 shows the computed attenuations of five body waves versus frequency and gas hydrate saturation at 50% porosity using the present model. For the cases of P1 and S1, as can be observed in Figure 5a,b, the attenuations of these two waves are minimal, both of them have attenuation peaks, and the attenuations of them initially increases, then decreases as frequency increases. Moreover, with the increase of hydrate saturation, the attenuation of P1 increases (Figure 5a), whereas that of S1 decreases at first but then increases (Figure 5b), which is consistent with the conclusion from the actual data [46]. For the cases of P2, S2 and P3, the three slow-wave attenuations are several orders of magnitude higher than attenuations of P1 and S1 (Figure 5c-e). Due to their slow propagation velocities and large attenuations in the seismic band, the three slow waves are difficult to be observed, which is similar to the slow-wave in Biot two-phase media. However, the three slow waves also carry part of the energy and decay quickly. Therefore, the effect of slow waves on attenuation should be considered when processing actual seismic data. With observations in Figures 4 and 5, we can expect that gas hydrate saturation can be estimated using our improved model.

Influence of Clay
The solid grain frame in porous media is usually composed of more than one substance, tiny amounts of clay have certain effects on wave propagation [52]. Clay usually exists in hydrate-bearing sediments, the deposition of gas hydrate comes under the influence of lithology greatly [52,53]. When clay exists in the solid grain frame, the bulk modulus and shear modulus can be expressed according to the Voigt-Reuss-Hill (VRH) average [54]: where K c and µ c define bulk modulus and shear modulus of clay, β represents the clay content in the solid grain frame. Generally, the clay content of a sand-rich reservoir with natural gas hydrate is small [34]. In this paper, a relationship between clay content and gas hydrate content is proposed, in which clay content is controlled within 20%: Figure 6 shows the velocities and attenuations calculated for a hydrate formation with the porosity φ = 0.5. For velocity, since the presence of clay makes the shear modulus of the frame decrease, the phase velocities of P1 and S1 are lower due to the occurrence of clay (Figure 6a,b). This phenomenon can be validated in the experiments by Han [55]. The velocity of P2 is not sensitive to clay content in the seismic band but sensitive in the acoustic logging band (Figure 6c), whereas the velocities of S2 and P3 are almost unchanged (Figure 6e,f). For attenuation, the absence of clay decreases the attenuation of P1 (Figure 6a), but doesn't change the attenuations of the other four waves (Figure 6b-e). Therefore, there is usually a certain amount of clay in the real hydrate reservoir, which has an impact on the phase velocity and attenuation that we should not neglect. The lack of clay will increase the saturation of inversion.

Model Validation
In this section, the effectiveness of the model can be verified by comparing the P1 velocity calculated by the improved Carcione-Leclaire model with the measured P-wave velocity of the Ocean Drilling Program (ODP) Leg 164 Site 955B. The leg found obvious BSR (Bottom Simulating Reflector) characteristic in the Black Ridge on the North American continental slope, and got hydrate samples. Note that, α is a function of sediment density, shear modulus, S-wave velocity, porosity, hydrate saturation and free parameter ε. ε has little influence on the calculation results and we take ε = 0.12, which was suggested by previous studies [46,47]. The S-wave velocity and porosity can be got from well log, while the hydrate saturation information can be obtained by pore water chlorinity measurement and methane content from pressure cores. Thus, the consolidation coefficient α of each location can be calculated using Equation (11).
The comparison between P1 velocity predicted by the theory in this paper and the actual acoustic well log data measured at the ODP Leg 164 Site 955B is shown in Figure 7.
The red dots represent the velocities calculated using the proposed model and the blue curve is the actual well log curve. The frequency is 20 kHz. The parameters used to calculate the consolidation coefficient are listed in Table 2. Figure 7 and Table 2 indicate that velocities calculated using the improved Carcione-Leclaire model fit well with the observed velocities. In particular, the maximum error between P1 velocities calculated by the model and P-wave velocities obtained by actual logging is 3.5%, and the average error of them is 1.6%. It indicates the theoretical model proposed in this paper has been proven to be effective. The difference could be partly due to the inaccurate sediment constituents used here. With observations in Figure 7, we can expect that P-wave velocity can be predicted according to this model for where it is missing. Conversely, the model can also be used to evaluate the logging data and estimate gas hydrate saturation. However, the densities and modulus used in our calculation are from the previous study [44], which will produce a certain error in the result. If we can get the composition of sediment and the content of each component accurately, the result will be more accurate and convincing.

Inversion of Hydrate Saturation
According to the analysis of ODP Leg 164 Site 955B, our proposed model can be used to quantitatively estimate the hydrate saturation of gas hydrate-bearing sediments using actual log data including the measured P-wave velocity and porosity information of ODP Leg 204 Site 1247B. The site is located on the west coast of Oregon in the Western Pacific Ocean of the United States, where there is a famous "Hydrate Ridge", which contains plenty of methane hydrate. Figure 8 shows the hydrate distribution calculated based on the improved Carcione-Leclaire model and acoustic logging data. The physical parameters used in the theoretical calculation are shown in Table 1, and the frequency selected in our calculation is 20 kHz. The most direct method is to deduce the calculation formula of hydrate saturation first, and then calculate the hydrate saturation by combining the measured velocity and porosity data. However, there is a certain degree of difficulty to solve the analytical solution of the gas hydrate saturation. Thus, the numerical solution is used to calculate the hydrate saturation in this section. First, we computed consolidation coefficients with various hydrate saturations from Equation (11), which used the measured S-wave velocities and porosities. Then the P-wave velocities are obtained by numerical calculation using our proposed model. When the calculated P-wave velocity is the same as the measured data, the hydrate saturation value used is the estimated saturation of the position. Figure 8 indicates that about 100-130m below the seafloor of Site 1247B is a hydrate-bearing layer, which is consistent with the analysis by Kuang [53]. The hydrate saturation of this layer is estimated as 0-27%, the upper limit is slightly higher than that estimated by Kuang. Further, the average saturation is calculated as 7.2%.

Conclusions
In this paper, an improved Carcione-Leclaire model is proposed to study wave propagation characteristics in gas hydrate-bearing sediments. The sensitivities of velocities and attenuations of P1 and S1 to various parameters especially gas hydrate saturation are discussed in detail. Further, the variation of phase velocities and attenuations of waves with frequency and gas hydrate saturation is emphatically studied. Moreover, the effects of clay on wave propagation characteristics are analyzed. Furthermore, we compared the velocity of P1 predicted by the model with the well log data from the Ocean Drilling Program (ODP) Leg 164 Site 995B to verify the applicability of the model. The findings are listed below: (1) Both velocities and attenuations of P1 and S1 are sensitive to gas hydrate saturation, especially attenuation. Furthermore, P1 and S1-velocity sensitivity to gas hydrate saturation increase when the basic value of gas hydrate saturation increases, and it is frequency independent. In the seismic band, P1 and S1-attenuation sensitivity decrease as the basic value of gas hydrate saturation increases while increases with the increase of frequency; In the ultrasonic band, P1 and S1-attenuation sensitivity increases with the increase of the basic value of gas hydrate saturation and decreases with the increase of frequency. Therefore, the smaller the gas hydrate saturation at low frequencies as well as the larger the gas hydrate saturation at high frequencies, the more reliable the inversion. (2) The numerical results show comprehensive relations of the velocity and attenuation with gas hydrate saturation and frequency. P1 and S1 are equivalent to P-and S-waves in single-phase elastic media which have almost no dispersion in the seismic band, while the velocities of them increase almost linearly as gas hydrate saturation increases. In the seismic band, P2, S2, and P3 are dispersive, and the attenuations are significant. In the ultrasonic band, however, P2, S2, and P3 are nondispersive and the attenuations are low. Besides, with the increase of gas hydrate saturation, the velocity of P2 initially decreases and then increases, the velocity of S2 increases, and the velocity of P3 initially increases and then decreases. Moreover, clay reduces the velocities of P1 and S1 significantly but has almost no effect on phase velocities of S2 and P3. The phase velocity of P2 is not sensitive to clay content in the seismic band, but sensitive in the ultrasonic band. The numerical simulations and field data analysis in this paper show that the wave velocity and attenuation characteristics are closely related to gas hydrate saturation, which can be used to estimate gas hydrate content. It also provides a potential for further applications in the evaluation of gas hydrate content via wave attenuation.
Author Contributions: The original idea was contributed by L.L. The manuscript preparation and simulation parts were finished by L.L. The follow-up guidance was given by X.Z. The supervisors of L.L. are X.W. and X.Z. The manuscript was improved and perfected by X.Z. and X.W. All authors have approved the publication of the paper. All authors have read and agreed to the published version of the manuscript.

Conflicts of Interest:
The authors declare no conflict of interest.

Abbreviations
List of Symbols b 12 friction coefficient between the solid grain frame and pore water b 13 friction coefficient between the solid grain frame and gas hydrate b 23 friction coefficient between the pore water and gas hydrate φ s volume fraction of solid grain φ w volume fraction of pore water φ h volume fraction of gas hydrate ρ s solid grain density ρ w water density ρ h gas hydrate density η w water viscosity κ s0 permeability of solid-grain frame κ h0 permeability of gas-hydrate frame K s solid grain bulk modulus K w water bulk modulus K h gas hydrate bulk modulus K hm bulk modulus of the gas hydrate frame K sm bulk modulus of the solid grain frame µ s solid grain shear modulus µ h gas hydrate shear modulus µ hm shear modulus of the gas hydrate frame µ sm shear modulus of the solid grain frame c 1 solid grain consolidation g 1 solid grain consolidation c 3 gas hydrate consolidation g 3 gas hydrate consolidation 21 tortuosity for water flowing through the solid grain frame a 23 tortuosity for water flowing through the gas hydrate a 23 = 1 + r 23 φ h (φ w ρ w + φ s ρ s )/φ w ρ w (φ w + φ s ) a 13 tortuosity for solid grain flowing through the gas hydrate a 31 tortuosity for gas hydrate flowing through the solid grain where m represents the number of constituents and f i denotes the volume fraction of component i. The bulk and shear moduli of sediment are calculated according to Gassmann equation [56]: The P-wave velocity can be expressed as: For non-contacting effective medium theory (EMT), gas hydrate floats in pores, it only changes the elastic properties of pore fluid but does not affect the frame properties. The effective bulk modulus of pore fluid is given by: (A5) For cementing EMT, hydrate is part of the frame, which reduces the original porosity, the porosity becomes φ = φ(1 − S h ), The modified bulk and shear moduli of effective minerals are given by: Besides, the natural gas hydrate is a component of minerals, the volume fraction of which can be expressed as: The three-phase Biot equation by Lee (TPBEL) is suitable for the calculation of velocity at low frequency, the P-wave velocity is given by [36]: where R i,j and ρ b are explained in Section 2.1.

Appendix B
EMT and TPBEL (see Appendix A) have been widely used to analyze the relationship between phase velocities and various parameters such as gas hydrate saturation [3,25,39,47]. In EMT, three microscopic modes of gas hydrate exist in hydrate-bearing sediments: hydrate is a pore fluid, which is called non-contacting EMT; gas hydrate cements sediment grains, which is named cementing EMT; hydrate contacts with sediment grains. A study has shown that non-contacting EMT and cementing EMT are major patterns of hydrate formation and decomposition [57]. Thus, we explore the variation of P-wave velocities (P1 in the improved Carcione-Leclaire model) with gas hydrate saturation in different models. As indicated in Figure A1, P-wave velocities calculated in TPBEL and improved Carcione-Leclaire model are obviously higher than two EMT patterns. Moreover, there is reasonable concordance between the P-wave velocities of TPBEL and the improved Carcione-Leclaire model at low frequency. However, at high frequency, when the hydrate saturation lower than 0.26 or higher than 0.71, the velocities predicted by the improved Carcione-Leclaire model are higher than those calculated by TPBEL. Figure A2 shows the variation of P-wave velocities with gas hydrate saturation in TPBEL and the improved Carcione-Leclaire model. The frequency has no effect on the velocities of TPBEL and EMT, but influence a little on the velocity of the improved Carcione-Leclaire model when the hydrate saturation lower than 0.26 or higher than 0.71. Thus, due to the TPBEL and EMT simplified the gas hydrate-bearing sediments, there is some error in the velocities of these models.