9 The mechanical behaviors of methane hydrate-bearing sands (MHBS) are largely a ected by the presence of methane hydrate, temperature, and pore pressure. In this study, we present a simple hypoplastic model for MHBS. Methane hydrate saturation is included as a state parameter a ecting the mechanical behaviors of MHBS. A new phase parameter is introduced to account for the coupled e ects of temperature and pore pressure on the me- chanical behaviors of MHBS. The phase parameter can be determined by a simple function of temperature and pore pressure. Comparison of the pre- dictions with experiments shows that the model is able to capture the salient behaviors of MHBS. Keywords: Hypoplastic constitutive model, Methane hydrate, Sands, 10 Phase state, Critical state 11 1. Introduction 12 Methane hydrate (MH) is a crystalline clathrate composed of methane 13 and water molecules. It is widely distributed beneath permafrost and in 14 marine continental margin sediments. MH can be used to extract methane 15 gas, and has attracted global interests as a promising energy resource [1]. 16 According to eld investigations, MHs primarily exist in the forms of discrete 17 nodules or veins in ne-grained sediments or disseminate in the pore spaces of 18 coarse sands. Compared to ne-grained sediments, sand reservoirs are more 19 Corresponding author. E-mail: shun.wang@boku.ac.at; mingjing.jiang@tju.edu.cnPreprint submitted to Acta Geotechnica November 13, 2019 attractive for MH exploitation due to the high-permeability of sand media 20 and the high-concentration of MHs [2, 3]. However, the exploitation of MH 21 may destabilize methane hydrate-bearing sands (MHBS) and potentially lead 22 to the unsteady deformation of seabed [4], in-stability of production well [5], 23 and even marine landslides [6]. It is therefore crucial to properly characterize 24 the mechanical behaviors of MHBS for MH exploitation. 25 The mechanical behaviors of MHBS have been extensively investigated by 26 laboratory tests. Sucient experimental results [7{13] show that the sti ness, 27 strength, and dilation of MHBS are enhanced by increasing MH saturation, 28 which is de ned as the ratio of MH volume to sediment pore volume. Besides, 29 the mechanical behaviors of MHBS are markedly in uenced by pore pressure 30 and temperature [7, 13]. For example, the peak strength of MHBS increases 31 by 18% when temperature decreases from 10C to 1C [7]; and a 40% increase 32 of the peak strength is observed when pore pressure increases from 5 MPa 33 to 10 MPa [13]. Hence, the MH saturation, temperature, and pore pressure 34 dependences are important in modeling the mechanical behaviors of MHBS. 35 During the last decade, a great e ort has been made in the development 36 of constitutive models for predicting the mechanical behaviors of MHBS. Ear- 37 lier models are usually developed based on the Ducan-Change (DC) theory 38 [14{16] and the Mohr-Coulomb (MC) theory [17{20]. In these constitutive 39 models, MH saturation and temperature are considered as important factors 40 that in uence the mechanical behaviors of MHBS. The DC type models es- 41 tablish a stress-strain relationship through a hyperbolic function, in which 42 the initial sti ness and the peak shear strength can be related to con ning 43 pressure [14], MH saturation [15], and temperature [14, 16]. With a nonlin- 44 ear elastic formulation, the nonlinear behavior of MHBS can be described, 45 while the applied hyperbolic function is usually too simple to fully capture 46 the volumetric and strain softening behaviors. In the MC type models, the 47 e ect of MH saturation can be considered by including the dependency of 48 cohesion on MH saturation [17{20]. However, without a yielding surface, the 49 MC type models predict inaccurate volumetric deformation of MHBS under 50 loading. 51 Although the aforementioned models are able to describe some basic fea- 52 tures of MHBS, the prediction of the critical state behavior, an important 53 aspect of the behavior of MHBS, is beyond the capacity of the models due to 54 the limitation of model framework. In contrast to these models, the critical 55 state models based on the Modi ed Cam-Clay framework consider the e ects 56 of MH saturation using the critical state parameters. The volumetric change, 57 2 strain softening, and MH-induced bonding can be reasonable described [21{ 58 25]. Thus, the critical state models are among the most common options in 59 predicting the mechanical behavior of MHBS. 60 A closer inspection of the relevant literature reveals, however, that the 61 e ects of temperature and pore pressure, have rarely been included in avail- 62 able critical state models. Although an attempt has been made to account 63 for the e ects of temperature and pore pressure within the framework of 64 bounding surface plasticity [23], the formulation of the constitutive equation 65 is complicated, which may hinder the understanding and application of the 66 model. From the practical point of view, a constitutive model with simple 67 formulation while is capable of capturing the salient behaviors of MHBS is 68 to be favoured. 69 Hypoplastic models describe the behaviors of granular materials by non- 70 linear tensorial functions [26]. Unlike elastoplastic models, hypoplastic mod- 71 els need not to prede ne yield surface, plastic potential, hardening rule, ow 72 rule, and the decomposition of deformation into elastic and plastic parts [26]. 73 Owing to their simple formulations with few model parameters, hypoplastic 74 models have been widely used to model the behaviors of sands [27{30], clay 75 [31, 32], frozen soil [33], and debris materials [34], etc. Recently, the hy- 76 poplstic model has also been applied to model the mechanical behaviors of 77 MHBS [35]. The e ects of MH saturation and con ning pressure on the 78 stress-strain and volumetric responses can be well captured. However, the 79 parameters for di erent MH saturations and con ning pressures are not uni- 80 ed, i.e., need to be recalibrated. Moreover, some salient factors, such as 81 void ratio, temperature, and pore pressure, are not discussed. 82 This study aims to develop a hypoplastic model accounting for the e ects 83 of MH saturation, stress level, void ratio, temperature, and pore pressure. 84 First, a simple hypoplastic model for sand [36] is enhanced by incorporating 85 the critical state concept. Second, the e ect of MH is incorporated by con- 86 sidering both the hydrate-induced bonding and the critical state line (CSL) 87 being dependent on MH saturation, temperature, and pore pressure. In par- 88 ticular, a phase parameter is introduced to account for the coupled e ects of 89 temperature and pore pressure on the mechanical behaviors of MHBS. After 90 that, the calibration of the model is illustrated. Last, the model is validated 91 by comparing simulation results with experiments. The performance of the 92 model and some limitations are also discussed. 93 3 2. Hypoplastic model with critical state for sand 94 Wu et al. [36] proposed a basic constitutive model for sand. This model 95 contains tensor polynomial with three linear terms and one nonlinear term 96 in strain rate: 97 =C1(tr)_+C2(tr)(tr_)I+C3tr(_) tr+C4(+)k_k (1) whereCi, (i= 1, 2, 3, 4) are dimensionless parameters, the sign tr denotes 98 the trace of a tensor, is the Cauchy stress tensor, _is the strain rate 99 (stretching) tensor, Iis the second-order unit tensor, is the deviatoric 100 stress tensor de ned as =(tr)I=3,k_kis the Euclidean norm of the 101 strain rate, and is the Jaumann stress rate de ned as: 102 =_+_!+_! (2) where _is the material time derivative of the Cauchy stress, and _!is the 103 spin tensor. 104 This basic model takes the stress tensor as the only state variable. As 105 a consequence, di erent parameters are required for the same material but 106 with di erent initial densities. This shortcoming, however, can be readily 107 remedied by integrating the critical state concept. To this end, a density 108 functionfd[37] is introduced to account for the e ect of void ratio. The 109 basic model (1) is then recast as follows: 110 =C1(tr)_+ (1fd)C2(tr)(tr_)I+C3tr(_) tr+fdC4(+)k_k(3) With the density function fd, this model is able to consider variation of 111 the density during shearing. Several forms of density functions can be found 112 in the literature. In this work, a simple form suggested by Wang et al. [38] 113 is used. It reads: 114 fd= (e ecrt) (4) whereeis the current void ratio, ecrtis the critical state void ratio, and is 115 a constitutive constant that controls the degree of strain softening [29, 38]. 116 The value of fdis less than 1 for a dense state, greater than 1 for a loose state, 117 and equals 1 at the critical state. The critical state void ratio is determined 118 4 by: 119 ecrt=eln(p0 pa) (5) whereis a parametric constants, p0= tr=3 is the e ective mean stress, pa 120 = 101.325 kPa is the atmospheric pressure, and eis the critical state void 121 at the atmospheric pressure. 122 The constitutive equation (3), however, does not account for the in uence 123 of the initial void ratio on the initial sti ness, and usually gives rise to a very 124 sti stress-strain curve in triaxial test simulations. To modulate the shape 125 of the stress-strain curve, the right side of Eq. (3) is multiplied by a sti ness 126 functionfsproposed by Wu and Kolymbas [36]: 127 fs=exp[ (ecrte)] (1 +r)2(6) where is a material parameter and ris the stress ratio kk=(tr). In this 128 way, the constitutive model for the host sand is formulated as: 129 =fs[C1(tr)_+(1fd)C2(tr)(tr_)I+C3tr(_) tr+fdC4(+)k_k] (7) 3. Model extension to MHBS 130 The mechanical behaviors of MHBS have been widely investigated by tri- 131 axial compression tests [7{13]. Experimental results reveal that the existence 132 of MH enhances the strength and dilation of MHBS. A higher MH saturation 133 usually leads to a larger cohesion and a higher location of the CSL on the 134 elogp0plane. In addition, the strength and dilatancy of MHBS can be 135 enhanced by a lower temperature and a higher initial pore pressure. 136 To account for the aforementioned e ects, we proceed to improve the 137 sand model (7) in two aspects. First, an internal variable Sis introduced 138 to account for MH-induced cohesion. Second, a new phase parameter L0is 139 proposed to describe the coupled e ects of temperature and pore pressure 140 on the cohesion and CSL of MHBS. The complete constitutive equations are 141 given in `Appendix' for clarity. 142 3.1. Incorporation of cohesion 143 Since the model (7) is homogeneous in stress, the limit surface in the 144 principal stress space is a cone with its apex at the origin. As a result, Eq. 145 5 (7) allows no tensile stress. This shortcoming can be removed by transforming 146 the stress tensor in Eq. (7) by an internal variable S: 147 T=+SI (8) whereTis a transformed stress tensor. In this way, the hydrate-induced 148 cohesion can be included by replacing byTin Eq. (7). 149 As indicated by experimental results [7, 8], MH saturation SMH, temper- 150 atureT, and pore pressure Phave signi cant in uences on the MH-induced 151 cohesion. Therefore, the value of Sshould be as a function of SMH,T, and 152 P, i.e.,S=f(SMH;T;P ). This function will be formulated in subsection 153 3.3.3 after introducing the new phase parameter L0. 154 3.2. The e ect of MH saturation 155 Experimental results reveal that the critical state line of MHBS will shift 156 upward on the elogp0plane with increasing MH saturation. Besides, the 157 CSL of MHBS with a given SMHis usually assumed to be parallel to that of 158 the host sand [7], as shown in Fig. 1. 159 /s48/s46/s49 /s49 /s49/s48 /s49/s48/s48/s48/s46/s53/s48/s46/s54/s48/s46/s55/s48/s46/s56/s48/s46/s57/s49/s46/s48 /s32/s32 /s73/s110/s99/s114/s101/s97/s115/s105/s110/s103/s32 /s83 /s77/s72/s86/s111/s105/s100/s32/s114/s97/s116/s105/s111/s32 /s101 /s77/s101/s97/s110/s32/s101/s102/s102/s101/s99/s116/s105/s118/s101/s32/s115/s116/s114/s101/s115/s115/s32 /s112/s39/s32/s72/s111/s115/s116/s32/s115/s97/s110/s100 /s32/s83 /s77/s72/s32/s61/s32/s50/s54/s46/s49/s37 /s32/s83 /s77/s72/s32/s61/s32/s51/s52/s46/s50/s37 /s32/s83 /s77/s72/s32/s61/s32/s52/s51/s46/s57/s37 /s32/s83 /s77/s72/s32/s61/s32/s53/s50/s46/s57/s37 /s32/s83 /s77/s72/s32/s61/s32/s54/s49/s46/s50/s37 Figure 1: Critical state lines for MHBS [7]. Based on the above experimental observations, the e ect of MH satura- 160 tion on the CSL is suggested as: 161 ecrt=e(MH)ln(p0=pa) (9) e(MH) =e+S 1 MH (10) 6 whereand 1are model parameters. The parameter eis the same as that 162 used in Eq. (5). 163 3.3. The e ects of temperature and pore pressure 164 Previous experimental works show that the changes of temperature and 165 pore pressure will lead to the variations of: (1) MH-induced cohesion and 166 (2) the location of the CSL of MHBS [7, 8]. In the following, we proceed 167 to consider the e ects of temperature and pore pressure in a uni ed way by 168 using a new phase parameter. 169 3.3.1. A new phase parameter L0170 The stability of MHs with respect to temperature and pore pressure can 171 be described by a phase parameter Lde ned by Hyodo et al. [39]. As illus- 172 trated in Fig. 2a, the phase parameter Lis de ned as the minimum distance 173 from a phase state ( ) to the stability boundary. This boundary di erenti- 174 ates the area where MHs are stable (solid) or unstable (gaseous). Generally, 175 MHs are stable only when the phase state lies above the stability boundary. 176 Within this area, a higher value of Lindicates that the MH is more stable, 177 and a higher strength can be gained [39], see Fig. 2b. Since the mechanical 178 behaviors of MHBS are predominately governed by the properties of MH, 179 the e ects of temperature and pore pressure on the mechanical behaviors of 180 MHBS can be considered using the phase parameter L[23, 40{43]. However, 181 the value of Lcan be obtained only by iterative method in previous studies 182 due to the high nonlinearity of the stability boundary. To facilitate the use 183 of the phase parameter concept, we will propose a simple stability boundary 184 based on experimental results [44{50]. 185 As shown in Fig. 3a, the stability boundary can be represented by two 186 lines in the T0=Tln(P=P 0) plane. This characteristic suggests that the new 187 phase parameter for a given phase state ( ) can be expressed as a piecewise 188 function as follows: 189 L0=8 >>>>>< >>>>>:j35(T0=T) + ln(P=P 0)35jp 352+ 12 2area ABD; j7:5(T0=T) + ln(P=P 0)7:5jp 7:52+ 12 2area CBD; 0 otherwise(11) 7 /s48/s46/s52 /s48/s46/s54 /s48/s46/s56 /s49/s46/s48 /s49/s46/s50/s48/s46/s48/s49/s48/s46/s49/s49/s49/s48/s49/s48/s48 /s32/s32/s79/s114/s105/s103/s105/s110/s97/s108/s32/s112/s104/s97/s115/s101/s32/s112/s97/s114/s97/s109/s101/s116/s101/s114/s32 /s76/s32/s32/s80/s104/s97/s115/s101/s32/s115/s116/s97/s116/s101/s32 /s32/s32/s40 /s80 /s44/s32 /s84 /s41/s32 /s32/s32 /s83/s116/s97/s98/s105/s108/s105/s116/s121/s32/s98/s111/s117/s110/s100/s97/s114/s121/s32/s111/s102/s32/s77/s72/s78/s111/s114/s109/s97/s108/s105/s122/s101/s100/s32/s112/s114/s101/s115/s115/s117/s114/s101/s32/s40 /s80 /s32/s91/s77/s80/s97/s93/s32/s47/s32/s50/s46/s54/s32/s77/s80/s97/s41 /s78/s111/s114/s109/s97/s108/s105/s122/s101/s100/s32/s116/s101/s109/s112/s101/s114/s97/s116/s117/s114/s101/s32/s40 /s84 /s32/s91/s75/s93/s32/s47/s32/s50/s55/s51/s32/s91/s75/s93/s41/s77/s72/s32/s100/s105/s115/s115/s111/s99/s105/s97/s116/s101/s100/s32/s40/s103/s97/s115/s101/s111/s117/s115/s41/s32/s97/s114/s101/s97(a) /s48/s46/s48/s48 /s48/s46/s48/s52 /s48/s46/s48/s56 /s48/s46/s49/s50 /s48/s46/s49/s54 /s48/s46/s50/s48/s48/s52/s56/s49/s50/s49/s54/s50/s48 /s32/s32 /s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s97/s108/s32/s100/s97/s116/s97 /s32/s70/s105/s116/s116/s105/s110/s103/s32/s108/s105/s110/s101/s77/s72/s32/s115/s116/s114/s101/s110/s103/s116/s104/s32/s40/s77/s80/s97/s41 /s65/s120/s105/s97/s108/s32/s115/s116/s114/s97/s105/s110/s32/s40/s37/s41 (b) Figure 2: (a) Normalized temperature and pressure stability boundary and (b) relationship between the strength of MHBS and L[39]. wherePandTdenote the current pore pressure and temperature, respec- 190 tively;P0= 2.23 MPa and T0= 273 K are the reference pore pressure and 191 temperature, respectively; area ABD and CBD are the acute angle areas 192 bounded by ABD and CBD, respectively (see Fig. 3b, in which line BD is 193 the angle bisector of ABC). Compared to the phase parameter given in pre- 194 vious literature [23, 40{43], the new phase parameter L0is much simpler and 195 can be explicitly determined. 196 /s48/s46/s56 /s49/s46/s48 /s49/s46/s50 /s49/s46/s52 /s49/s46/s54/s45/s54/s45/s52/s45/s50/s48/s50/s52/s54 /s32/s68/s117/s97/s110/s32/s101/s116/s32/s97/s108/s46/s32/s91/s52/s52/s93/s32/s32/s32/s32/s32/s32 /s32/s77/s97/s107/s111/s103/s111/s110/s32/s97/s110/s100/s32/s83/s108/s111/s97/s110/s32/s91/s52/s53/s93 /s32/s74/s97/s103/s101/s114/s32/s97/s110/s100/s32/s83/s108/s111/s97/s110/s32/s91/s52/s54/s93 /s32/s77/s97/s101/s107/s97/s119/s97/s32/s91/s52/s55/s93 /s32/s83/s117/s110/s32/s101/s116/s32/s97/s108/s46/s32/s91/s52/s56/s93/s32 /s32/s83/s97/s98/s105/s108/s32/s101/s116/s32/s97/s108/s46/s32/s91/s52/s57/s93 /s32/s67/s104/s97/s32/s101/s116/s32/s97/s108/s46/s32/s91/s53/s48/s93 /s32/s32 /s32/s70/s105/s116/s116/s105/s110/s103/s32/s108/s105/s110/s101 /s108/s110/s40 /s80 /s47/s80 /s48/s41/s32/s61/s32/s45/s55/s46/s53/s215 /s40 /s84 /s48/s47/s84 /s41/s43/s55/s46/s53/s108/s110/s40 /s80 /s47/s80 /s48/s41/s32 /s84 /s48/s47/s84 /s32/s108/s110/s40 /s80 /s47/s80 /s48/s41/s32/s61/s32/s45/s51/s53/s215 /s40 /s84 /s48/s47/s84 /s41/s43/s51/s53 (a) /s48/s46/s57/s48 /s48/s46/s57/s53 /s49/s46/s48/s48 /s49/s46/s48/s53 /s49/s46/s49/s48 /s49/s46/s49/s53/s45/s49/s48/s49/s50/s51 /s76 /s39/s32/s32/s32/s32/s32/s32/s32/s32/s32/s32/s32 /s32/s32/s40 /s80 /s44/s32 /s84 /s41/s32/s68 /s67/s66/s65 /s32/s32/s108/s110/s40 /s80 /s47/s80 /s48/s41 /s84 /s48/s47/s84 (b) Figure 3: Stability boundary of methane hydrate in pure water: (a) experimental data and (b) the determining of the new phase parameter L0. 8 3.3.2. The e ect of L0on the CSL 197 The relationship between the CSL of MHBS and the phase parameter L0198 can be established by assuming that the phase parameter changes only the 199 position of the CSL rather than its slope (Fig. 4). The parameter controlling 200 the position of the CSL of MHBS is considered to be dependent on the phase 201 parameterL0. The CSL of MHBS thus has the following formulation: 202 ecrt=e(L0)ln(p0=pa) (12) e(L0)=e+L0 2(13) where 2is an additional parameter controlling the e ect of L0on the CSL 203 of MHBS. 204 /s48/s46/s49 /s49 /s49/s48 /s49/s48/s48/s48/s46/s53/s48/s46/s54/s48/s46/s55/s48/s46/s56/s48/s46/s57/s49/s46/s48 /s83 /s77/s72/s32/s61/s32/s52/s51/s46/s57/s37 /s73/s110/s99/s114/s101/s97/s115/s105/s110/s103/s32 /s76/s39 /s32/s32 /s32/s72/s111/s115/s116/s32/s115/s97/s110/s100 /s32/s76/s39 /s32/s61/s32/s48/s46/s48/s50/s52/s57 /s32/s76/s39 /s32/s61/s32/s48/s46/s48/s48/s53/s49/s86/s111/s105/s100/s32/s114/s97/s116/s105/s111/s32 /s101 /s77/s101/s97/s110/s32/s101/s102/s102/s101/s99/s116/s105/s118/s101/s32/s115/s116/s114/s101/s115/s115/s32 /s112/s39 /s32/s40/s77/s80/s97/s41 Figure 4: The in uence of phase parameter L0on the CSL of MHBS (data from Hyodo et al. [7]) 3.3.3. The e ects of SMHandL0on the CSL and cohesion 205 Based on Eqs. (10) and (13), we proposed the following functions to 206 account for the e ects of MH saturation and phase state on the CSL of 207 MHBS: 208 ecrt=e(MH;L0)ln(p0=pa) (14) e(MH;L0)=e+S 1 MHL0 2(15) 9 In addition, the e ects of MH saturation and phase state on the MH- 209 induced cohesion is suggested as: 210 S=aSb MHL0(16) whereaandbare additional parameters controlling the e ects of the new 211 phase parameter and MH saturation on the MH-induced cohesion. 212 4. Calibration of the model 213 The proposed model for MHBS contains 14 parameters, 8 of which de- 214 scribe the mechanical behaviors of the host sand and the remaining 6 account 215 for the e ects of MH saturation and the phase state. The value of the phase 216 state parameter L0can be explicitly obtained by Eq. (11) for desired tem- 217 perature and pore pressure. For instance, Table 1 shows the values of L0for 218 desired temperature and pore pressure. Other parameters are determined in 219 the following subsections. 220 Table 1: The values of L0for desired temperature and pore pressure Temperature T[K] Pore pressure P[MPa] State parameter L0 274 10 0.0392 274 15 0.0508 278 5 0.0051 278 10 0.0249 283 10 0.0075 4.1. Calibration of parameters for the host sand 221 The parameters Ci, (i=1, 2, 3, 4) can be determined with a single triaxial 222 compression test according to the procedure described by Wu and Bauer 223 [27]. In this study, the material parameters determining Ci, (i=1, 2, 3, 224 4) areEi/c= 180,i= 0.2, crt= 0, and crt= 31[7], whereEiis 225 the initial tangent modulus; c= 1 MPa is the con ning pressure in the 226 triaxial compression test; iis the initial Poisson ratio; crtis the critical 227 state dilatancy angle; and crtis the critical state friction angle. The critical 228 state parameters eandcan be determined by tting the critical state line 229 of the host sand [7], as shown in Fig. 4. 230 10 The parameters and can be calibrated by tting the peak strength 231 and modulus of MHBS (in the same single triaxial compression test). The 232 parameters for modeling of the host sand [7] are shown in Table 2. 233 Table 2: Model parameters for the host sand Parameter Value C1 -50 C2 -53.46 C3 -436.48 C4 -147.73 e 1.145  0.126 0.8 3 * and should be re-calibrated for MHBS. 4.2. Calibration of parameters for MHBS 234 Apart from the 8 model parameters for the host sand, 6 additional pa- 235 rameters are used to account for the e ects of MH saturation, temperature, 236 and pore pressure. As shown in Eqs. (15) and (16), both the MH saturation 237 and the phase state in uence the hydrate-induced cohesion and CSL on the 238 elnp0plane. For the MH-induced cohesion, the internal variable Scan be 239 determined by back calculating the experimental data of the peak strength 240 of MHBS [21]. As illustrated in Fig. 5, the internal variable Sis taken as the 241 mean e ective stress p0corresponding to zero deviatoric stress q(q=0 10 3, 242 with0 1and0 3being the maximum and minimum principal stresses, respec- 243 tively) in the peak strength envelope. Therefore, the parameters aandbcan 244 be obtained by tting the variation of SagainstSMHandL0. 245 As illustrated in Fig. 6, the value of e(MH;L0)in Eq. (15) is the intercept 246 of the CSL on the y-axis. The parameters , 1, and 2can be obtained by 247 tting the variation of e(MH;L0)againstSMHandL0. Note that the parame- 248 ters and need to be re-calibrated for MHBS, which may be because the 249 existence of MH alters the particle-scale property of sand grains. By tting 250 a single triaxial compression test on MHBS, the values of and are deter- 251 mined as 0.3 and 5, respectively. The remaining 5 parameters are presented 252 in Table 3. 253 11 /s45/s50 /s48 /s50 /s52 /s54 /s56 /s49/s48 /s49/s50/s48/s51/s54/s57/s49/s50/s49/s53 /s68/s101/s99/s114/s101/s97/s115/s105/s110/s103/s32 /s83 /s77/s72/s32/s97/s110/s100/s32 /s76/s39 /s32/s32 /s32/s83 /s77/s72/s32/s61/s32/s53/s51/s46/s55/s37/s44/s32 /s76/s39 /s32/s61/s32/s48/s46/s48/s50/s52/s57/s32 /s32/s83 /s77/s72/s32/s61/s32/s53/s51/s46/s50/s37/s44/s32 /s76/s39 /s32/s61/s32/s48/s46/s48/s53/s48/s56/s32 /s32/s83 /s77/s72/s32/s61/s32/s52/s48/s46/s51/s37/s44/s32 /s76/s39 /s32/s61/s32/s48/s46/s48/s50/s52/s57 /s32/s83 /s77/s72/s32/s61/s32/s52/s50/s46/s57/s37/s44/s32 /s76/s39 /s32/s61/s32/s48/s46/s48/s48/s53/s49 /s32/s32/s32/s32/s32/s32/s32/s32/s32/s32/s70/s105/s116/s116/s105/s110/s103/s32/s108/s105/s110/s101/s68/s101/s118/s105/s97/s116/s111/s114/s32/s115/s116/s114/s101/s115/s115/s32 /s113/s32 /s40/s77/s80/s97/s41 /s69/s102/s102/s101/s99/s116/s105/s118/s101/s32/s109/s101/s97/s110/s32/s115/s116/s114/s101/s115/s115/s32 /s112/s39 /s32/s40/s77/s80/s97/s41Figure 5: Peak strength envelopes of MHBS at di erent hydrate saturations (data from Hyodo et al. [7]). Table 3: The remaining parameters for MHBS Parameter a[MPa]b  1 2 Value 359 5.815 6.526 2.486 0.523 5. Model performance 254 In this section, model predictions are compared with experimental results 255 of MHBS. A series of numerical triaxial compression tests, considering dif- 256 ferent con ning pressures, initial void ratios, MH saturations, temperatures, 257 and pore pressures, are carried out to validate this model. The parameters 258 used in the simulations are summarized in Tables 2 and 3. 259 5.1. Model performance for the host sand 260 Figures 7 and 8 show the comparisons between the predictions and ex- 261 periments on host sand with di erent con ning pressures and void ratios, 262 respectively. In experiments, the strength and sti ness of MHBS increase 263 with the increase of the e ective con ning pressure, while the dilatancy is 264 weakened by the increase of the void ratio. In addition, a strain-softening 265 behavior is observed under a lower e ective con ning pressure, i.e., 1 MPa, 266 while strain-hardening behaviors are observed under higher e ective con- 267 ning pressures of 3 and 5 MPa. Overall, the model can well predict the 268 12 /s48/s46/s49 /s49 /s49/s48 /s49/s48/s48/s48/s46/s53/s48/s46/s54/s48/s46/s55/s48/s46/s56/s48/s46/s57 /s32/s32 /s32/s32/s32/s32/s32/s32/s32/s32/s32/s32/s70/s105/s116/s116/s105/s110/s103/s32/s108/s105/s110/s101 /s32/s72/s111/s115/s116/s32/s115/s97/s110/s100 /s32/s83 /s77/s72/s32/s61/s32/s50/s53/s46/s54/s44/s32 /s76/s39 /s32/s61/s32/s48/s46/s48/s50/s52/s57 /s32/s83 /s77/s72/s32/s61/s32/s51/s54/s46/s57/s44/s32 /s76/s39 /s32/s61/s32/s48/s46/s48/s50/s52/s57 /s32/s83 /s77/s72/s32/s61/s32/s52/s51/s46/s53/s44/s32 /s76/s39 /s32/s61/s32/s48/s46/s48/s50/s52/s57 /s32/s83 /s77/s72/s32/s61/s32/s53/s51/s46/s51/s44/s32 /s76/s39 /s32/s61/s32/s48/s46/s48/s50/s52/s57 /s32/s83 /s77/s72/s32/s61/s32/s52/s52/s46/s51/s44/s32 /s76/s39 /s32/s61/s32/s48/s46/s48/s48/s53/s49/s73/s110/s99/s114/s101/s97/s115/s105/s110/s103/s32 /s83 /s77/s72/s32/s97/s110/s100/s32 /s76/s39/s86/s111/s105/s100/s32/s114/s97/s116/s105/s111/s32 /s101 /s77/s101/s97/s110/s32/s101/s102/s102/s101/s99/s116/s105/s118/s101/s32/s115/s116/s114/s101/s115/s115/s32 /s112/s39 /s32/s40/s77/s80/s97/s41Figure 6: The e ects of SMHandL0on the critical state line (data from Hyodo et al.[7]). stress-strain and volumetric responses of the host sand for di erent con ning 269 pressures and void ratios observed in experiments [7]. 270 /s48 /s51 /s54 /s57 /s49/s50 /s49/s53/s48/s51/s54/s57/s49/s50/s49/s53/s68/s101/s118/s105/s97/s116/s111/s114/s105/s99/s32/s115/s116/s114/s101/s115/s115/s32 /s113 /s32/s40/s77/s80/s97/s41/s32/s101 /s48/s32/s61/s32/s48/s46/s54/s54 /s32/s32 /s65/s120/s105/s97/s108/s32/s115/s116/s114/s97/s105/s110/s32/s40/s37/s41/s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44/s32 /s99/s39/s32/s61/s32/s49/s32/s77/s80/s97 /s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44/s32 /s99/s39/s32/s61/s32/s51/s32/s77/s80/s97 /s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44/s32 /s99/s39/s32/s61/s32/s53/s32/s77/s80/s97 /s32/s77/s111/s100/s101/s108/s32/s112/s114/s101/s100/s105/s99/s116/s105/s111/s110 (a) /s48 /s51 /s54 /s57 /s49/s50 /s49/s53/s53/s48/s45/s53 /s32/s101 /s48/s32/s61/s32/s48/s46/s54/s54 /s32/s32 /s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44/s32 /s99/s39/s32/s61/s32/s49/s32/s77/s80/s97 /s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44/s32 /s99/s39/s32/s61/s32/s51/s32/s77/s80/s97 /s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44/s32 /s99/s39/s32/s61/s32/s53/s32/s77/s80/s97 /s32/s77/s111/s100/s101/s108/s32/s112/s114/s101/s100/s105/s99/s116/s105/s111/s110/s86/s111/s108/s117/s109/s101/s116/s114/s105/s99/s32/s115/s116/s114/s97/s105/s110/s32/s40/s37/s41 /s65/s120/s105/s97/s108/s32/s115/s116/s114/s97/s105/s110/s32/s40/s37/s41 (b) Figure 7: Comparisons between predicted and experimental results for the host sand with di erent e ective con ning pressures: (a) stress-strain and (b) volumetric responses (data from Hyodo et al. [7]). 5.2. Model performance for MHBS 271 Figure 9 presents the comparisons between experimental and the pre- 272 dicted results of MHBS with di erent MH saturations. Experimental results 273 13 /s48 /s51 /s54 /s57 /s49/s50 /s49/s53/s48/s50/s52/s54/s56/s49/s48 /s99/s39/s32/s61/s32/s51/s32/s77/s80/s97 /s32/s32/s68/s101/s118/s105/s97/s116/s111/s114/s105/s99/s32/s115/s116/s114/s101/s115/s115/s32 /s113 /s32/s40/s77/s80/s97/s41 /s65/s120/s105/s97/s108/s32/s115/s116/s114/s97/s105/s110/s32/s40/s37/s41/s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44/s32 /s101 /s48/s32/s61/s32/s48/s46/s54/s54 /s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44/s32 /s101 /s48/s32/s61/s32/s48/s46/s56/s49 /s32/s77/s111/s100/s101/s108/s32/s112/s114/s101/s100/s105/s99/s116/s105/s111/s110(a) /s48 /s51 /s54 /s57 /s49/s50 /s49/s53/s56/s54/s52/s50/s48/s45/s50 /s32/s32 /s99/s39/s32/s61/s32/s51/s32/s77/s80/s97/s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44/s32 /s101 /s48/s32/s61/s32/s48/s46/s54/s54 /s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44/s32 /s101 /s48/s32/s61/s32/s48/s46/s56/s49 /s32/s77/s111/s100/s101/s108/s32/s112/s114/s101/s100/s105/s99/s116/s105/s111/s110/s86/s111/s108/s117/s109/s101/s116/s114/s105/s99/s32/s115/s116/s114/s97/s105/s110/s32/s40/s37/s41 /s65/s120/s105/s97/s108/s32/s115/s116/s114/s97/s105/s110/s32/s40/s37/s41 (b) Figure 8: Comparisons between predicted and experimental results for the host sand with di erent void ratios: (a) stress-strain and (b) volumetric responses (data from Hyodo et al. [7]). indicate that, with the increase of SMH, the MHBS samples show enhance- 274 ment in strength, sti ness, and dilatancy. At an MH saturation of 53.1%, a 275 signi cant dilative behavior is observed regardless of the adopted high e ec- 276 tive con ning pressure (5 MPa). Despite the slightly discrepancy, the model 277 successfully captures the main e ects of SMHon the stress-strain behaviors 278 and the volumetric responses of MHBS. 279 /s48 /s51 /s54 /s57 /s49/s50 /s49/s53/s48/s51/s54/s57/s49/s50/s49/s53/s68/s101/s118/s105/s97/s116/s111/s114/s105/s99/s32/s115/s116/s114/s101/s115/s115/s32 /s113 /s32/s40/s77/s80/s97/s41 /s32/s32 /s99/s39/s32/s61/s32/s53/s32/s77/s80/s97 /s76/s39 /s32/s61/s32/s48/s46/s48/s50/s52/s57/s32 /s65/s120/s105/s97/s108/s32/s115/s116/s114/s97/s105/s110/s32/s40/s37/s41/s32/s72/s111/s115/s116/s32/s115/s97/s110/s100 /s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44/s32 /s83 /s77/s72/s32/s61/s32/s50/s52/s46/s50/s37 /s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44/s32 /s83 /s77/s72/s32/s61/s32/s51/s53/s46/s49/s37 /s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44/s32 /s83 /s77/s72/s32/s61/s32/s53/s51/s46/s49/s37 /s32/s77/s111/s100/s101/s108/s32/s112/s114/s101/s100/s105/s99/s116/s105/s111/s110 (a) /s48 /s51 /s54 /s57 /s49/s50 /s49/s53/s53/s48/s45/s53 /s32/s32 /s99/s39/s32/s61/s32/s53/s32/s77/s80/s97 /s76/s39/s32 /s61/s32/s48/s46/s48/s50/s52/s57/s32/s86/s111/s108/s117/s109/s101/s116/s114/s105/s99/s32/s115/s116/s114/s97/s105/s110/s32/s40/s37/s41 /s65/s120/s105/s97/s108/s32/s115/s116/s114/s97/s105/s110/s32/s40/s37/s41/s32/s72/s111/s115/s116/s32/s115/s97/s110/s100/s32/s32/s32/s32/s32/s32/s32/s32/s32 /s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44/s32 /s83 /s77/s72/s32/s61/s32/s50/s52/s46/s50/s37 /s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44/s32 /s83 /s77/s72/s32/s61/s32/s51/s53/s46/s49/s37 /s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44/s32 /s83 /s77/s72/s32/s61/s32/s53/s51/s46/s49/s37 /s32/s77/s111/s100/s101/s108/s32/s112/s114/s101/s100/s105/s99/s116/s105/s111/s110 (b) Figure 9: Comparisons between predicted and experimental results for MHBS with di erent MH saturations: (a) stress-strain and (b) volumetric responses (data from Hyodo et al. [7]). Figures 10 and 11 illustrate the predicted mechanical behaviors of MHBS 280 14 with di erent con ning pressures and void ratios, respectively. It can be ob- 281 served in experiments that for MHBS with SMHof 53.7%, the strength, sti - 282 ness, and dilatancy are enhanced by increasing e ective con ning pressures 283 and decreasing void ratios. The proposed model can capture the changes 284 of the stress-strain behaviors and volume change responses with varying 285 e ective con ning pressures and void ratios. Admittedly, there are some 286 discrepancies between the predicted stress-strain/volumetric responses and 287 experimental data. The discrepancies maybe attributed to the e ect of MH 288 morphology, which is not considered in the current model. This e ect will 289 be discussed in Section 6. 290 /s48 /s51 /s54 /s57 /s49/s50 /s49/s53/s48/s53/s49/s48/s49/s53/s50/s48 /s32/s32 /s83 /s77/s72/s61/s32/s53/s51/s46/s55/s37 /s76/s39 /s32/s61/s32/s48/s46/s48/s50/s52/s57/s32/s68/s101/s118/s105/s97/s116/s111/s114/s105/s99/s32/s115/s116/s114/s101/s115/s115/s32/s40/s77/s80/s97/s41 /s65/s120/s105/s97/s108/s32/s115/s116/s114/s97/s105/s110/s32/s40/s37/s41/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44 /s32 /s99/s39/s32/s61/s32/s49/s32/s77/s80/s97 /s32/s32 /s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44 /s32 /s99/s39/s32/s61/s32/s51/s32/s77/s80/s97 /s32/s32 /s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44 /s32 /s99/s39/s32/s61/s32/s53/s32/s77/s80/s97 /s32 /s32/s77/s111/s100/s101/s108/s32/s112/s114/s101/s100/s105/s99/s116/s105/s111/s110 (a) /s48 /s51 /s54 /s57 /s49/s50 /s49/s53/s52/s50/s48/s45/s50/s45/s52/s45/s54/s45/s56/s45/s49/s48 /s83 /s77/s72/s61/s32/s53/s51/s46/s55/s37 /s76/s39 /s32/s61/s32/s48/s46/s48/s50/s52/s57/s32 /s32/s32/s86/s111/s108/s117/s109/s101/s116/s114/s105/s99/s32/s115/s116/s114/s97/s105/s110/s32/s40/s37/s41 /s65/s120/s105/s97/s108/s32/s115/s116/s114/s97/s105/s110/s32/s40/s37/s41/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44 /s32 /s99/s39/s32/s61/s32/s49/s32/s77/s80/s97 /s32/s32 /s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44 /s32 /s99/s39/s32/s61/s32/s51/s32/s77/s80/s97 /s32/s32 /s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44 /s32 /s99/s39/s32/s61/s32/s53/s32/s77/s80/s97 /s32 /s32/s77/s111/s100/s101/s108/s32/s112/s114/s101/s100/s105/s99/s116/s105/s111/s110 (b) Figure 10: Comparisons between predicted and experimental results for MHBS with di erent con ning pressures: (a) stress-strain and (b) volumetric responses (data from Hyodo et al. [7]). Figures 12 and 13 compare the measured and predicted results of MHBS 291 under di erent temperatures and pore pressures. Experimental results show 292 that when pore pressure increase from 5 MPa to 10 MPa, the sti ness of 293 MHBS increase markedly; and the strength and dilatancy are also enhanced. 294 In addition, when temperature decreases from 283 K to 274 K, obvious in- 295 crease in the sti ness, strength, and dilatancy of MHBS are observed. The 296 comparisons of experimental and numerical results show that the proposed 297 model is able to suciently capture the above e ects. Overall, the results 298 predicted by our model are close to the experimental data. 299 15 /s48 /s51 /s54 /s57 /s49/s50 /s49/s53/s48/s51/s54/s57/s49/s50/s49/s53 /s83 /s77/s72/s61/s32/s53/s50/s37 /s76/s39 /s32/s61/s32/s48/s46/s48/s50/s52/s57/s32 /s32/s32/s68/s101/s118/s105/s97/s116/s111/s114/s32/s115/s116/s114/s101/s115/s115/s32/s40/s77/s80/s97/s41 /s65/s120/s105/s97/s108/s32/s115/s116/s114/s97/s105/s110/s32/s40/s37/s41/s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44/s32 /s101 /s48/s32/s61/s32/s48/s46/s56/s49 /s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44/s32 /s101 /s48/s32/s61/s32/s48/s46/s54/s54 /s32/s77/s111/s100/s101/s108/s32/s112/s114/s101/s100/s105/s99/s116/s105/s111/s110(a) /s48 /s51 /s54 /s57 /s49/s50 /s49/s53/s53/s48/s45/s53 /s83 /s77/s72/s61/s32/s53/s50/s37 /s76/s39 /s32/s61/s32/s48/s46/s48/s50/s52/s57/s32 /s32/s32/s86/s111/s108/s117/s109/s101/s116/s114/s105/s99/s32/s115/s116/s114/s97/s105/s110/s32/s40/s37/s41 /s65/s120/s105/s97/s108/s32/s115/s116/s114/s97/s105/s110/s32/s40/s37/s41/s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44/s32 /s101 /s48/s32/s61/s32/s48/s46/s56/s49 /s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44/s32 /s101 /s48/s32/s61/s32/s48/s46/s54/s54 /s32/s77/s111/s100/s101/s108/s32/s112/s114/s101/s100/s105/s99/s116/s105/s111/s110 (b) Figure 11: Comparisons between predicted and experimental results for MHBS with di erent void ratios: (a) stress-strain and (b) volumetric responses (data from Hyodo et al. [7]). /s48 /s51 /s54 /s57 /s49/s50 /s49/s53/s48/s50/s52/s54/s56/s49/s48/s49/s50 /s32/s32 /s83 /s77/s72/s61/s32/s52/s48/s46/s56/s37/s68/s101/s118/s105/s97/s116/s111/s114/s105/s99/s32/s115/s116/s114/s101/s115/s115/s32/s40/s77/s80/s97/s41 /s65/s120/s105/s97/s108/s32/s115/s116/s114/s97/s105/s110/s32/s40/s37/s41/s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44/s32/s112/s111/s114/s101/s32/s112/s114/s101/s115/s115/s117/s114/s101/s32/s111/s102/s32/s53/s32/s77/s80/s97 /s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44/s32/s112/s111/s114/s101/s32/s112/s114/s101/s115/s115/s117/s114/s101/s32/s111/s102/s32/s49/s48/s32/s77/s80/s97 /s32/s77/s111/s100/s101/s108/s32/s112/s114/s101/s100/s105/s99/s116/s105/s111/s110 (a) /s48 /s51 /s54 /s57 /s49/s50 /s49/s53/s53/s48/s45/s53 /s32/s32/s86/s111/s108/s117/s109/s101/s116/s114/s105/s99/s32/s115/s116/s114/s97/s105/s110/s32/s40/s37/s41 /s65/s120/s105/s97/s108/s32/s115/s116/s114/s97/s105/s110/s32/s40/s37/s41/s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44/s32/s112/s111/s114/s101/s32/s112/s114/s101/s115/s115/s117/s114/s101/s32/s111/s102/s32/s53/s32/s77/s80/s97 /s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44/s32/s112/s111/s114/s101/s32/s112/s114/s101/s115/s115/s117/s114/s101/s32/s111/s102/s32/s49/s48/s32/s77/s80/s97 /s32/s77/s111/s100/s101/s108/s32/s112/s114/s101/s100/s105/s99/s116/s105/s111/s110/s83 /s77/s72/s61/s32/s52/s48/s46/s56/s37 (b) Figure 12: Comparisons between predicted and experimental results for MHBS with di erent pore pressures: (a) stress-strain and (b) volumetric responses (data from Hyodo et al. [7]). 6. Discussion 300 The mechanical behaviors of MHBS show a strong dependency on MH 301 saturation. With the increase of MH saturation, a signi cant change in me- 302 chanical response will occur. Besides, the coupled e ects of temperature 303 and pore pressure are also important aspects of the mechanical behaviors of 304 MHBS. So far, most of available constitutive models for MHBS take only the 305 e ect of MH saturation into account. One exception is the contribution made 306 16 /s48 /s51 /s54 /s57 /s49/s50 /s49/s53/s48/s51/s54/s57/s49/s50/s49/s53 /s32/s32 /s83 /s77/s72/s32/s61/s32/s53/s50/s46/s52/s37/s68/s101/s118/s105/s97/s116/s111/s114/s105/s99/s32/s115/s116/s114/s101/s115/s115/s32/s40/s77/s80/s97/s41 /s65/s120/s105/s97/s108/s32/s115/s116/s114/s97/s105/s110/s32/s40/s37/s41/s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44/s32 /s84 /s32/s61/s32/s50/s56/s51/s32/s75 /s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44/s32 /s84 /s32/s61/s32/s50/s55/s52/s32/s75 /s32/s77/s111/s100/s101/s108/s32/s112/s114/s101/s100/s105/s99/s116/s105/s111/s110(a) /s48 /s51 /s54 /s57 /s49/s50 /s49/s53/s53/s48/s45/s53 /s32/s32/s86/s111/s108/s117/s109/s101/s116/s114/s105/s99/s32/s115/s116/s114/s97/s105/s110/s32/s40/s37/s41 /s65/s120/s105/s97/s108/s32/s115/s116/s114/s97/s105/s110/s32/s40/s37/s41/s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44/s32 /s84 /s32/s61/s32/s50/s56/s51/s32/s75 /s32/s69/s120/s112/s101/s114/s105/s109/s101/s110/s116/s115/s44/s32 /s84 /s32/s61/s32/s50/s55/s52/s32/s75 /s32/s77/s111/s100/s101/s108/s32/s112/s114/s101/s100/s105/s99/s116/s105/s111/s110/s83 /s77/s72/s32/s61/s32/s53/s50/s46/s52/s37 (b) Figure 13: Comparisons between predicted and experimental results for MHBS with di erent temperatures: (a) stress-strain and (b) volumetric responses (data from Hyodo et al. [7]). by Ng et al. [23], who consider the coupled e ects of temperature and pore 307 pressure by a phase parameter similar to that in this study. However, the 308 phase parameter in [23] is determined in an indirect way (e.g., by iterative 309 method), while the phase parameter proposed in this study can be deter- 310 mined straightforward by a simple function. In addition, the formulation of 311 the proposed model is much more simple than that of most available consti- 312 tutive models for MHBS. Therefore, the proposed hypoplastic model possess 313 great advantages in modeling the salient mechanical behaviors of MHBS. 314 Besides, the proposed phase parameter may provide convenience in future 315 modeling of the coupled e ects of temperature and pore pressure. 316 Apart from MH saturation, temperature, and pore pressure, the MH mor- 317 phology is also an important factor in controlling the mechanical behaviors 318 of MHBS. The MH morphology is usually categorized into four types: pore- 319 lling, load-bearing, cementation, and grain-coating [12]. In the past, the 320 MH morphology is viewed to be controlled by hydrate formation environ- 321 ment. For example, in a free water environment, hydrate tends to form from 322 dissolved gas molecules or gas bubbles, resulting in the pore- lling type or 323 the load-bearing type; in a free gas environment, hydrate tends to form as the 324 cementation type. Based on this idea, experiments performed in di erent hy- 325 drate formation environments (free gas/free water environment) [7, 8] were 326 usually used for the validation of constitutive models in terms of di erent 327 MH morphologies (cementation/pore- lling type) [23, 24, 51]. However, MH 328 17 morphology cannot be identi ed simply by the hydrate formation environ- 329 ment. Based on high-quality tomography (e.g., micro-computed tomography 330 and X-ray tomography), it is reported that pore- lling type MH is formed 331 even in the free gas environment [52]. Moreover, reliable evidences that the 332 transformation of MH morphology and coexistence of di erent MH morpholo- 333 gies are increasingly reported [53, 54]. Therefore, cautions should be taken 334 when using experimental data in terms of MH morphology since idealized 335 pore- lling or cementation type MH may not physically exist. 336 The e ect of MH morphology is not accounted for in this study because 337 available experiments are not sucient to isolate the e ects of pore- lling, 338 cementation, grain-coating, and load-bearing types of MH. Inevitably, the 339 ignorance of MH morphology hinders the accuracy in model prediction. In 340 addition, the accuracy in model prediction is also in uenced by the adopted 341 critical state parameter. Since the critical state parameters are determined 342 at the end of shearing, where the critical state may not be reached [7], the 343 performance of our model may be a ected. To establish a comprehensive 344 constitutive model for MHBS, the e ect of MH morphology and the critical 345 state behaviors of MHBS need further study. 346 7. Conclusions 347 A state-dependent critical state model is developed within the frame- 348 work of hypoplasticity. The model is able to account for hydrate saturation, 349 temperature, and pore pressure. A new phase parameter with explicit for- 350 mulation is proposed, which includes the coupled e ects of temperature and 351 pore pressure on the mechanical behaviors of MHBS in a uni ed way. The 352 model is validated with experimental results of MHBS obtained under a wide 353 range of stress level, void ratios, hydrate saturation, temperature, and pore 354 pressure. Comparisons of the predictions with experiments show the model 355 is satisfactory in predicting the salient mechanical behaviors of MHBS. 356 Bene ted from the non-linear tensorial functions, the hypoplastic model 357 has the advantage of simplicity while retaining its ability to capture the 358 salient behaviors of MHBS. These advantages make the hypoplastic model 359 as an attractive alternative to the prevailing elastoplastic constitutive models 360 for MHBS. 361 18 Acknowledgments 362 The authors would like to acknowledge the National Natural Science 363 Foundation of China with Grant Numbers 51639008 and 51890911, the H2020 364 Marie Sk lodowska-Curie Actions RISE 2017 \HERCULES" with grant num- 365 ber 778360, and \FRAMED" with grant number 734485. The second author 366 wishes to thank the Otto Pregl Foundation for nancial support in Austria. 367 Appendix 368 The complete hypoplastic equations for MHBS are illustrated as follows: =fsh C1(trT)_+(1fd)C2(trT)(tr_)I+C3tr(T_) trTT+fdC4(T+ T)k_ki where 8 >>>>>>< >>>>>>:T=+SI S=aSb MHL0 fs=exp[ (ecrte)] (1+r)2 fd= (e ecrt) ecrt=e(MH;L0)ln(p0=pa) e(MH;L0)=e+S 1 MHL0 2 In the above equations, Tis the transformed stress tensor, is the 369 Jaumann stress rate, is the Cauchy stress tensor, _is the stretching tensor, 370 Iis the second-order unit tensor, Sis the internal variable eis the current void 371 ratio,ecrtis the critical state void ratio, ris the stress ratio kk=(tr),L0is 372 the phase state parameter, e(MH;L0)is the critical state void of MHBS at the 373 atmospheric pressure, eis the critical state void of sand at the atmospheric 374 pressurepa, andp0is the e ective mean stress. The parameters of the model 375 areC1,C2,C3,C4,e,, , ,a,b,, 1, 2, andL0. 376