赵亚鹏 刘乐乐 孔 亮†,3) 刘昌岭 吴能友

* (青岛理工大学土木工程学院,山东青岛 266033)

† (中国地质调查局青岛海洋地质研究所,自然资源部天然气水合物重点实验室,山东青岛 266071)

** (青岛海洋科学与技术国家实验室海洋矿产资源评价与探测技术功能实验室,山东青岛 266071)

†† (青岛理工大学理学院,山东青岛 266033)

引 言

天然气水合物(简称水合物)是由天然气分子与水分子在高压、低温环境下形成的非化学计量笼型结晶物质[1],广泛分布于全球近海大陆架及陆地永久冻土地带[2],资源储量巨大,被认为是21 世纪最具潜力的替代能源之一[3].自20 世纪开始,世界各国陆续成立了专门的水合物研究机构并开展了相关研究[4],在此基础上,前苏联、加拿大、美国及日本等国家[5-8]相继进行了水合物试采,证明了水合物开采的技术可行性.与此同时,经过20 余年的不懈努力,我国在南海海域发现了储量高达800 亿吨油当量的水合物资源[9],并于2017 和2020年在南海北部海域两次成功试采水合物[10-11],创造了多项世界记录,率先实现了泥质粉砂型水合物资源的安全可控开采.

水合物资源的开发前景主要受开采技术条件及开发成本的控制,应以实际工程地质条件为基础,结合储层特性选取合理的开采方式[12].当前水合物开采方法主要有热激、降压、注入化学试剂、二氧化碳置换以及多方法联合等[13-17],虽然开采方法各异,但是均需原位分解水合物[18-19].水合物分解导致沉积物颗粒胶结程度降低,并可能产生超静孔压,储层力学强度随之弱化,进而引发各种工程地质灾害,如海底滑坡、海底浊流、井壁失稳及平台倾覆等[20-23].因此,厘清开采过程中含水合物土的力学特性演化规律是实现工程地质风险可控的重要理论基础,是实现水合物资源安全、高效、持续开采的前提条件之一[24-25].

国内外学者围绕含水合物土的力学特性开展了大量的研究工作,研究态势呈现出多尺度融合、多方法联用、多学科交叉等特征.目前实验研究主要以高压低温水合物三轴剪切实验为主[26],辅以电学、光学、声学等测量手段[27-29],再结合理论分析、数值模拟[30-31]等研究方法,系统研究含水合物土的力学特性演化机制.研究表明,含水合物土的力学特性由其内部结构控制,特别是水合物的微观行为扮演了极其重要的角色,越来越多的学者着眼于水合物微观力学领域[32],试图从微观结构演化解释宏观力学现象,“宏微观相结合”的科学范式成为当前乃至未来含水合物土力学特性研究的主流[33].因此,本文从水合物晶体(分子尺度)、水合物与土颗粒界面(颗粒尺度)、含水合物土(岩心尺度)3 个尺度对含水合物土微观力学特性研究现状进行综述,归纳微观力学特性研究的最新进展,阐述主要的微观测试技术,梳理出当前含水合物土微观力学研究的不足并给出针对性建议,以期为含水合物土力学特性研究提供新的思路与方向.

1 水合物晶体力学性质研究

1.1 晶体结构类型

水合物晶体的“主体”水分子之间由氢键连接并形成“笼型”结构,而“客体”天然气分子则通过范德华力填充于“笼”中[34].学者们现已发现多达100 余种水合物,按结构将其分为3 类,即sI,sII 及sH 型[35].不同类型水合物具有不同的晶体结构,在抗压强度和体积模量等力学性质上也表出明显的差异[36-37].sI,sII 和sH 型水合物所能容纳的烃类气体分子依次增大,而最终形成何种类型水合物则与气体组分及其含量有关.不同结构类型水合物对比可见表1.

表1 不同类型天然气水合物晶体结构对比[36,38]Table 1 Comparison of different structures in natural gas hydrate crystals[36,38]

1.2 研究方法

1.2.1 实验研究

考虑到晶粒尺寸对多晶固体力学性能的显著影响[39],学者们预猜想多晶水合物在力学特性上也将受到晶粒大小的影响[40].针对该问题,早期研究对比分析了海洋水合物与合成水合物的晶粒尺寸分布[41],并进行了毫米尺度下二氧化碳(CO2)水合物、甲烷(CH4)水合物以及四氢呋喃(C4H8O)水合物的抗拉强度研究[42],印证了晶粒尺寸确实会对水合物晶体力学特性尤其是塑性应变产生重要的影响.Atig 等[43]介绍了一种非常规的水合物微观力学特性研究方法,探讨了多晶甲烷水合物在微米尺度下的拉伸特性.该研究首先在薄玻璃毛细管中生成水合物,然后将非接触的热应力施加于水合物外壳(图1),进而研究了微米级别的多晶水合物抗拉特性.研究结果表明,温度和退火时间可控制水合物晶粒大小,水合物脆性随退火时间的延长而增强,而极限强度及应变则显著降低;颗粒结构的冰晶体并不是多晶甲烷水合物的理想相似物,海域甲烷水合物的抗拉强度可能远低于目前业界的认识.

图1 非接触应力施加步骤[43]Fig.1 Principal application steps of contactless stress[43]

水合物颗粒间的黏聚力与周围水的接触时间、接触力及界面能成正比,并随温度的升高而增大[44],水合物颗粒之间的液桥力是引起颗粒黏附力的主要原因[45].当含水合物土受到外力作用而发生颗粒旋转、滚动时,首先需要克服颗粒液桥力,这一点在水合物力学特性的原子力显微镜(atomic force microscopy,AFM)研究中已得到证明.液桥力示意图及AFM 实验液桥模型如图2 所示.

图2 液桥力示意图及AFM 实验液桥模型[46-47]Fig.2 Schematic of liquid bridge force and the model of liquid bridge for AFM experiment[46-47]

根据相关文献,即使在低于相平衡的温度下,水合物表面仍然存在纳米级别的似液层[48]以及压力作用下的分解液[49],气体水合物与砂颗粒之间亦会存在[50].颗粒间不同作用力(范德华力、液桥力、流体作用力等)的理论计算[51],则从定量的角度上进一步证明了液桥力是导致水合物颗粒黏附力的主要原因.而颗粒间液桥力受多因素的影响,以温度及液桥体积影响为切入点,对环戊烷水合物颗粒间的相互作用研究表明[52],未转化的水滴对颗粒团簇具有重要影响,液桥的形成能显著增强颗粒间黏附力,而黏附力、接触面积则与液滴体积没有直接的关系;同时,水合物颗粒呈现出与冰颗粒相似的力学行为[53].借助于水合物AFM,学者们用液桥力的公式进行了计算[46],发现探针与水合物间黏附力随温度的降低而减小.此外,受环境因素的影响,水合物颗粒与液体水滴呈动态作用关系,两者不断转化的同时伴随体积的持续变化,这将对颗粒间的液桥力产生显著影响[54].Aman 等[45,47,55]开展了较为系统的水合物液桥力研究,通过改进的微机械力学装置提出了一种直接测量颗粒间黏附力的新方法,对气−水两相黏附力进行对比分析,进一步探讨了不同添加剂对水合物颗粒之间液桥力的影响.

现阶段颗粒液桥研究主要集中于原子力显微镜等微观装置的实验分析,在进一步扩大研究尺度以及宏观联系方面范围有限,应逐步开展颗粒液桥的宏观力学影响研究.另一方面,虽然目前已经建立起了一定数量的液桥力微观理论模型,但是考虑多因素影响的理论模型仍相对匮乏,也有待进一步地深入研究.

1.2.2 分子动力学模拟

在晶体力学性质研究方面,除实验研究外,分子动力学(molecular dynamics,MD)模拟也是较为常用的研究方法.MD 以牛顿力学为基础,结合叠加原理,能够对分子尺度的微观力学进行模拟研究.

MD 研究表明水合物晶体强度、弹性模量及变形等与客体分子的大小、形状、极性等密切相关[56],同时还受到应变速率、温度、客体分子占有率的影响[57].水合物晶体与冰晶体的力学特性相近,一般表现为脆性破坏,而水合物断裂取决于客体分子类型,变形则受氢键的角位移控制.主客体分子关系的MD 研究发现[58],甲烷分子在“大笼”及“小笼”中均呈现各向同性,二氧化碳分子则在“大笼”中表现出明显的各向异性;而在混合气体中,无论是客体分子之间,还是主客体分子之间,甲烷的存在都会导致二氧化碳分子分布倾向性的改变.与纯水合物相比,冰含量对于含冰水合物强度影响较大,水合物−水合物双晶结构抗剪强度远大于水合物−冰双晶抗剪强度[59].在单晶及多晶水合物力学特性方面,单晶水合物呈现脆性破坏,多晶则为延性破坏,多晶甲烷水合物的力学稳定性与颗粒大小及形态密切相关,这种依赖关系归因于晶粒尺寸对晶界变形的影响[60].此外,学者们围绕天然气水合物晶体,通过MD 方法开展了包括不同气体分子置换过程[61]、水合物生长成核[62-63]、分解动力学[64]及分子作用模型[65]等方面的研究,进一步补充和丰富了水合物晶体分子尺度研究.水合物的典型MD 模拟如图3 所示.

图3 甲烷水合物分子结构[59]Fig.3 Molecular structure of methane hydrate[59]

相对而言,MD 模拟是一种非常规的研究手段,可以从分子及原子尺度进行微观、介观到宏观的物理力学性质解释.在水合物微观力学特性研究尤其是在水合物成核、分解及热力学方面,起到了极大地推动作用.而当前水合物MD 模拟研究主要基于均匀成核,而真实自然环境中存在非均匀成核情况,如何合理模拟非均匀成核带来的影响将是MD 研究的一大挑战.近期水合物晶体MD 模拟研究汇总可见表2.

表2 近期水合物晶体MD 模拟研究汇总Table 2 Summary of recent MD research on hydrate crystals

2 水合物与土颗粒界面力学性质研究

将水合物晶体结构尺度进一步扩大,晶体团簇形成颗粒,水合物颗粒与沉积物颗粒之间按照接触形态的不同可分为4 种主要赋存模式[67](如图4 所示),图4 中GH 为水合物,SP 为土颗粒.图4(a)为孔隙填充模式.水合物在沉积物颗粒孔隙之间生成,与沉积物颗粒无实质接触[68],此时水合物饱和度较低;图4(b)为颗粒支撑模式.在颗粒填充模式基础上,水合物颗粒进一步扩大并与沉积物颗粒接触,共同承担应力传递,此时水合物饱和度达到25%～40%[69];图4(c)为胶结模式.与颗粒填充模式不同,胶结模式水合物首先在沉积物颗粒接触处生成,并将周围沉积物颗粒“黏结”在一起,形成整体骨架[70];图4(d)为包裹模式.该模式下水合物将沉积物颗粒完全包裹.

图4 粗粒土中水合物赋存模式[67]Fig.4 Idealized gas hydrate morphologies in coarse-grained soils[67]

以上4 种模式主要存在于粗粒沉积物中,除孔隙填充模式外,其余3 种模式均会对力学特性产生显著影响[71].在“限气”或“富水”情况下易形成填充模式,在“限水”或“富气”情况下易形成包裹及胶结模式[72-73].细粒沉积物中则以结核状、脉状、块状、透镜状等模式存在[74-75](见图5),生成过程倾向于“挤开”沉积物发生“颗粒替代”,这与粗粒沉积物区别明显.上述赋存模式可进一步分为粗粒沉积物中均匀分布的孔隙填充型水合物以及细粒沉积物中非均匀分布的裂隙填充型水合物[76-77].对于含不同赋存类型水合物土的力学性质,很大程度上受到水合物与土颗粒界面力学行为的影响与控制.

图5 含水合物细粒土实物岩心[77]Fig.5 Hydrate-bearing fine-grained soils in nature[77]

2.1 研究方法

2.1.1 模拟实验研究

水合物与土颗粒界面力学性质的实验研究以计算机断层扫描(computed tomography,CT)联用的三轴剪切实验为代表,还包括部分基于微观实验装置、微观测试技术的颗粒界面关系研究.Lei 等[78]的CT 联用三轴实验表明,水合物主要通过限制颗粒旋转强化颗粒骨架,进而实现沉积物强度及刚度的增大,而应变软化的主要原因在于水合物与砂颗粒之间裂纹的产生;剪切应力会引起颗粒破碎并产生破坏裂缝,水合物则“流入”剪切裂缝或被尖锐的破碎颗粒所刺破.针对剪切过程的局部变形问题[79-80],借助于CT 图像,学者们确定了试样离散位置的局部刚度,并发现局部刚度与水合物饱和度成正相关关系,其中不含水合物的沉积物样品表现出延性破坏,含水合物样品则呈现出脆性破坏.基于微米尺度图像的局部变形(剪切面、颗粒变形及旋转等)量化分析[81],发现剪切带内沉积物的结构变化是基质颗粒和水合物颗粒旋转运动的结果,剪切带厚度随水合物饱和度的增大而减小.此外,水合物首先在颗粒表面发生分解,随后向颗粒孔隙发生转移;而水合物分解会引起沉积物组构的变化,进一步增强砂土颗粒分布的各向异性,并造成颗粒支撑及胶结效应的显著降低,颗粒间的均匀分布则可显著增强砂粒体系的力链稳定性[82].表3 给出了近期水合物CT 三轴实验汇总.

表3 近期含水合物土CT 联用三轴剪切实验汇总Table 3 Summary of current triaxial shearing tests combined with CT conducted on hydrate-bearing soils

除CT 三轴实验外,基于微观力学测量装置的直接测量结果表明[44,86],界面黏附力与颗粒接触时间及载荷大小有关,在接触时间一定的条件下,黏附力随载荷的增大而增大;颗粒摩擦特性方面,水合物颗粒与沉积物颗粒之间摩擦系数在不同载荷下呈现出不同的变化规律,当载荷较小时,随载荷增大,水合物颗粒与沉积物颗粒间摩擦系数亦不断增大.根据水合物的AFM 研究[87],颗粒间摩擦特性的差异在于颗粒表面性质或颗粒接触关系的不同,如颗粒表面粗糙度、颗粒尺寸、颗粒间似液层厚度等;水合物的不定形特征主要受低温及小尺寸晶粒的诱导,水合物表面形貌及粗糙度受环境温度、接触介质属性及晶粒尺寸所控制.不同个体分子的水合物抗拉实验则进一步表明[88],水合物在基质上成核以及水合物与基质整体接触时,表面粗糙程度对宏观力学特性的影响是不同的.此外,与水合物晶体颗粒类似,水合物与沉积物颗粒间黏聚力也会受到液桥力的影响,且沉积物基质表面对液桥力具有积极影响[45].

2.1.2 数值模拟研究

与有限元解决连续问题不同,离散元法(discrete element method,DEM)的颗粒流程序(particle flow code,PFC)不受变形量限制,能够对剪切过程中的开裂、分离等现象进行模拟,特别适用于解决非连续介质力学问题[89],被广泛应用于微观力学分析,并被认为是联系宏−微观力学的有效桥梁[32].

根据PFC 研究,颗粒界面间的黏结作用是含水合物土宏观力学响应的微观本质,这种黏结作用会对受力过程中的颗粒滚动产生显著影响,限制颗粒翻滚,进而影响宏观性质.CT 与扫描电子显微镜(scanning electron microscope,SEM)研究[90-91]表明,沉积物颗粒之间具有不同的颗粒形状,呈现不规则特性,颗粒形状对土体的影响则主要体现在法向接触力偏离形心引起的滚动力矩作用[92].同时不同形状颗粒具有不同的抗滚动系数,基于颗粒抗滚动作用的微观数值模拟表明[93],抗滚动系数的增大虽然会导致配位数的降低、孔隙率的增大,但是由于平均接触力的增大,颗粒之间将会形成更强的力链,进而引起宏观力学特性的增强.目前室内实验大多采用气饱和法合成水合物,生成的水合物以颗粒包裹型为主,相关的PFC 数值模拟研究表明,从颗粒界面关系的角度出发,包裹型水合物存在3 种微观作用机制[94],即颗粒黏结、颗粒扩大以及颗粒成角(黏结断裂时,存在粗糙断裂面,增大颗粒粗糙度).基于3 种微观机制的真三轴DEM 模拟[95]表明,包裹型水合物颗粒间的黏结作用相对较弱,真三轴剪切过程的力学特性接近于无黏结颗粒松散材料.此外,水合物剪切过程中在剪切带内外存在明显的颗粒演化差异[96],剪切带内部黏结大量破坏,颗粒发生转动,局部孔隙率不断增大;剪切带外部黏结破坏则相对较少,颗粒几乎不发生转动.这也进一步表明水合物与沉积物颗粒界面间的黏结关系是水合物存在而导致的微观效应,其会造成沉积物剪切破断过程中的微观颗粒旋转抑制,并呈现区域差异性,最终导致沉积物宏观力学特性的演变.

PFC 对水合物力学特性的有效模拟依赖于其内部的计算模型,而水合物颗粒与基质颗粒之间的界面作用关系,尤其是黏结效应则是PFC 计算模型建立的基础.不同类型水合物与沉积物颗粒之间具有不同的黏结厚度,这种界面关系将会呈现出宏观上的差异性特征[97].针对该问题,借助于SEM 图片,可获得水合物饱和度与黏结厚度之间的函数关系,进而建立考虑水合物黏结厚度的微观力学胶结模型[98].在多场耦合研究方面,基于实验结果,还可建立考虑应力、温度影响的黏结接触离散元模型[99],能够有效反映应变局部化特征;以及考虑水合物黏结厚度的热与流体力学耦合模型[100],对该模型的双轴压缩实验表明,随着颗粒间黏结厚度的增大,水合物峰值强度及剪胀率均显著增大.将计算流体动力学与DEM 相结合,学者们提出了一种考虑水合物速率依赖性的黏结接触模型[101](如图6 所示),基于该模型的不排水循环剪切实验表明,力链、接触组构以及平均旋转速率等微观变量与宏观力学特性之间具有紧密的联系.此外,部分学者认为水合物的存在会导致摩擦性的增强,而不一定与黏结有关,据此提出了含水合物土的非黏结微观模型[102].模型中将水合物表示为固体颗粒,即使在小应变下也会产生骨架效应.颗粒粒径与颗粒摩擦力之间存在着平衡关系,利用与孔径分布相关的附加因子可以实现两者的统一.

2.2 研究结果

基于实验结果,学者们试图从水合物与沉积物颗粒间界面关系,尤其是其在剪切过程中的动态演化规律对宏观力学特性进行解释.现有实验数据表明,水合物的存在会对沉积物硬化、软化、剪缩、剪胀等特性产生显著影响,颗粒尺度分析认为这是水合物颗粒、沉积物颗粒及孔隙空间相对关系的宏观体现[103-104].

水合物的生成会占据孔隙空间,这种“填充”作用导致沉积物整体骨架致密性的提高,同时水合物对砂颗粒具有显著的“黏结”作用,两者共同对力学特性产生影响.而不同赋存模式下含水合物土的“填充”与“黏结”作用也会存在明显差异,其受水合物饱和度的控制,内在本质则是微观组构关系.针对剪切过程中的含水合物土,随着剪切应力的升高,颗粒间接触增多,并不断摩擦、挤压及旋转,最终发生颗粒破碎以及“黏结”破坏.即剪切过程的实质是“填充”与“黏结”作用控制下的颗粒关系重新排列.对于密实沉积物而言,往往发生剪胀,松散沉积物则对应剪缩[105].

剪切过程的颗粒尺度作用机理可见图7.当水合物饱和度较低时,水合物对沉积物的强度影响较小,此时水合物与砂颗粒之间接触较少,“黏结”作用不明显[106-107],而水合物自身的“胶结”强度很容易被剪切应力所突破[108].随着水合物饱和度的增大,水合物与砂颗粒发生接触,“黏结”作用增强,水合物成为颗粒骨架并与砂颗粒共同承担应力传递.此时水合物的存在会显著限制基质颗粒在剪切应力下的滚动及重新排列.在此基础上,水合物强度演化受水合物“胶结”强度、砂颗粒与水合物之间“黏结”强度相对大小的控制.当“黏结”强度大于“胶结”强度时,首先在水合物中产生剪切面,沉积物强度演化由水合物力学特性控制,这种情况很可能发生在表面粗糙或具有较低比表面积沉积物中;当“胶结”强度大于“黏结”强度时,水合物与砂颗粒界面处将产生剪切面,并在剪切面形成大量水合物颗粒,此种情况很可能发生在表面光滑或具有较高比表面积沉积物中[69].

图7 颗粒尺度下含水合物土剪切机理[105]Fig.7 Particle-scale shearing mechanism of hydrate-bearing soils[105]

针对水合物颗粒的剥离与破碎问题,一些学者认为含水合物土在拉伸及剪切过程中分别具有不同的破坏模式[109].一种假想的微观机制如图8 所示,当受拉应力作用时,颗粒剥离或拉伸破坏主要取决于界面性质,并呈现出两种破坏模式,分别为沉积物颗粒与水合物颗粒之间的剥离以及水合物自身内部的断裂.失效的水合物与固体颗粒之间没有进一步的相互作用,两种模式如图8(a)所示.而在剪切破坏中,破碎颗粒将继续与其他颗粒产生相互作用,宏观上表现为峰后强度[69,110-111].与拉伸破坏相比,剪切破坏存在混合机制,即颗粒剥离过程中伴随水合物的断裂,相比纯颗粒剥离而言其体变更小.

图8 含水合物土拉伸及剪切破坏机制(改自文献[109])Fig.8 Failure modes of tension and shearing in hydrate-bearing soils (modified from Ref.[109])

近期,Wu 等[85]通过CT 三轴实验对水合物的拉伸及剪切过程进行了更深一步的探讨分析.CT 图像(图9)表明,无论是哪种类型的水合物,在剪切过程中,水合物颗粒均呈现整体移动,颗粒间相对位置关系几乎不变,而不含水合物的砂颗粒即使在较小轴向应变下也呈现无序移动.这无疑进一步印证了水合物对沉积物整体强度的显著增强作用.并且水合物对颗粒重新排列的抑制作用在应力应变的线性区域(εa≤2%)效果明显;随着轴向应变增大,应力应变曲线进入塑性及屈服阶段,水合物脱落、颗粒破碎,并出现了明显的剪切带,最终在剪切作用下发生了颗粒的重新排列.

图9 含水合物土剪切过程中截面CT 图像[85]Fig.9 CT images of hydrate-bearing soils during shearing[85]

3 含水合物土力学性质研究

3.1 离散元数值模拟

PFC 对水合物的有效模拟以制样方法为基础,重点在于突出水合物及沉积物颗粒之间的相互作用关系,进而获取这种微观关系所导致的宏观力学特性.就黏结效应而言,现有研究发现微观黏结参数对宏观力学特性的影响存在差异[112],随黏结强度的增大,峰值强度及弹性模量均增大;随黏结半径的增大,弹性模量增大而峰值强度基本不变.在抗滚动方面,由于黏结存在而导致的抗滚动作用会显著提高含水合物土的抗剪强度以及剪胀性,但是其存在一个阈值,当超过阈值后,这种增强作用将逐渐减小[93];并且不同的抗滚动作用会形成不同的微观力链分布(如图10 所示),这也将对沉积物宏观力学特性产生显著影响.含水合物土的PFC2D 数值模拟表明[113],剪切带的形成与黏结的大量破坏有关,并与应变局部化、微观变量局部化密切相关.以黏土质沉积物为研究对象,并考虑颗粒间毛细水效应的模拟研究发现[114],含水合物土具有剪缩及应变硬化的力学特性;随水合物饱和度的增大,颗粒间毛细水合物将产生一种键合效应,其随颗粒间隙的增大而线性减小.根据多类型实验(三轴压缩、各向同性压缩、恒应力比压缩及真三轴压缩)的水合物宏微观力学特性研究[115],颗粒包裹型水合物力学特性表现为弱黏结砂,颗粒胶结在加载初期就破裂,而在加载后基质砂颗粒的摩擦特性随颗粒形貌的变化而变化.对循环载荷下的含水合物土DEM 模拟表明[116],沉积物的破坏模式取决于循环载荷幅值以及水合物饱和度,组构的各向异性存在阈值,超过阈值沉积物结构则会趋于破坏.此外,针对含水合物土的黏聚力及内摩擦角力学参数,PFC 模拟研究发现[117],黏聚力和内摩擦角随剪胀量的增大而减小,应变局部化与黏结破坏、接触力链、颗粒速度等微观局部变量密切相关.当水合物饱和度增大时,水合物抗剪强度、黏聚力及体缩均明显增大,而摩擦角基本不变[96,100,118],这与宏观三轴实验结果一致[119].

图10 不同抗滚动系数下含水合物土内力链分布[93]Fig.10 Force chain distribution in hydrate-bearing soils under different rolling resistance coefficients[93]

微观参数对PFC 模拟结果至关重要,但是由于离散元颗粒与真实颗粒之间的差异,微观参数的确定较为困难,且不能直接获得.目前主要通过“试算”的方式确定微观参数[89,112],首先以宏观水合物三轴实验为依据,在一定范围内赋予接触模型微观力学参数,根据模拟结果不断调试,最终获得有效反应宏观力学特性的微观力学参数[120],据此开展一系列的微观数值模拟研究.同时基于前述的PFC 内部计算模型分析,可知,构建有效反应水合物力学特性的内部计算模型是水合物PFC 研究的另一难点与基础.针对该问题,蒋明镜等[121-122]开展了早期的探索研究,并进行了有关模型建立与参数反演的深入分析,研究成果对水合物PFC 研究起到了很好的推动作用.

DEM 以PFC 为代表在岩石力学及岩土工程领域得到了广泛应用,当前水合物力学研究所涉及的诸多力学问题都可划归为岩土力学问题[33].同时PFC 可实现水合物从微米颗粒尺度[123]到海底滑坡千米尺度[124]的有效模拟,因此可以预见,未来在水合物宏−微观相结合研究尤其是微观力学特性方面,PFC 还将继续发挥出更大的作用.

3.2 本构模型研究

本构模型研究主要以力学实验结果为依据,结合弹塑性本构建模理论及相关假设进行模型建立.在水合物方面,目前较多的是基于水合物存在对沉积物宏观力学特性的影响而开展相关研究,例如通过含水合物土的三轴剪切实验,对不同水合物饱和度的应力应变关系进行分析,进而建立起适用于中细粒含水合物土的本构模型[125].微观力学本构则主要围绕水合物与沉积物之间的微观作用关系、力学作用机理等进行模型构建,如考虑水合物赋存模式、填充及黏结效应、基质颗粒作用关系、颗粒破碎及损伤等.

水合物的微观赋存模式对含水合物土的宏观力学特性具有显著影响,这已经成为普遍共识.基于这一考虑角度,通过引入水合物“有效饱和度”的概念,并重新定义有效应力,可以建立起考虑水合物微观赋存模式的弹塑性本构模型[126],能够有效反映水合物微观赋存模式及水合物饱和度对含水合物土宏观力学特性的影响,对刚度、强度、应变软化等力学特性的演变均能较好地描述[127].针对水合物的微观黏结机理,在经典的临界状态模型框架内,通过引入非线性黏结与线性脱结规律,并将临界状态线表示为水合物饱和度的函数,建立起考虑水合物微观黏结机理以及黏结强度的本构模型[128],能够充分反映水合物饱和度、围压以及密度对含水合物土力学特性的影响;此外,引入次加载概念,用增强因子表征水合物的微观胶结效应,同时结合颗粒损伤,则可建立起能够有效预测不同孔隙习性(如孔隙填充、孔隙胶结等)水合物力学行为的新型本构模型[129].水合物的存在会造成含水合物土微观孔隙结构的改变,部分学者将这种孔隙改变归结为填充及胶结效应,通过引入压硬性参量,建立起考虑水合物填充及胶结效应的弹塑性本构模型[130].沉积物基质与水合物颗粒之间存在明显的界面关系,这将造成宏观受力过程中,微观力学响应的差异,即两者具有不同的力学贡献,据此,将宏观受力进行微观划分,同时引入“分区应力”的概念,则可建立起能够有效描述水合物加载及分解过程中的本构模型[131].Fang 等[132]根据塑性滑移理论,将宏观受力分解为体积响应以及空间分布的一系列微剪切响应(如图11 所示),进而得到了本构方程,该模型可以综合描述含水合物土的固结、硬化、软化、剪胀及非共轴等特性.

图11 含水合物土的多剪切模型概念框架[132]Fig.11 Illustration of conceptual framework for multishear model of hydrate-bearing soils[132]

针对水合物开采出砂问题,通过不同方法对微观颗粒运移的本质进行量化研究,发现水合物的存在减缓了应力松弛幅度,据此提出了考虑应力松弛特性的本构模型[133],该模型可实现水合物储层产气行为的有效模拟.部分学者[134]认为含水合物土缺乏真正的黏聚力,其力学行为受孔隙尺度的颗粒运动所控制.据此建立起了一种全新的水合物弹塑性本构模型,该模型并不依赖于水合物晶体与沉积物颗粒间的物理黏结,而是依赖于水合物侵入孔隙并对沉积物产生的致密化效应.根据该模型,水合物形成过程中沉积物微观孔隙体积的减小使其结构变硬,并随着沉积物密度的增大而产生近似的力学效应.此外一些学者们从损伤力学的角度出发,认为应力−应变关系的塑性阶段是内部颗粒发生损伤并逐渐累积的过程[135],损伤过程则可分为微观结构的快速损伤与完全损伤阶段[136],同时结合微元强度、细观力学机制等微观假设,建立起了一系列的基于损伤理论的含水合物土本构模型[137-139].其建模思路及建模过程对水合物微观力学本构模型的建立起到了很好的借鉴作用.

水合物本构模型及微观力学本构模拟方面已经积累了大量的研究成果,而其主要问题在于模型适用性及参数确定方面,常用的模型参数确定方法包括实验获取、拟合分析以及合理假设等[127,130,140-142].绝大多数模型基于特定实验条件或特定假设而建立,适用范围有限.在模型参数方面,部分模型的参数众多,逐一确定参数取值十分困难,很多模型参数并不具有实际物理意义,这无疑进一步增大了确定参数的难度.从建模角度出发,上述本构模型可进一步划分为3 类,分别为弹塑性模型、统计损伤模型以及其他模型.而宏观水合物本构模型还包括非线性弹性模型,然而由于该类模型采用弹性框架,在水合物微观力学特性建模方面受到了较大局限,报道较少.含水合物土微观本构模型分类及主要特点如表4 所示.

表4 含水合物土微观本构模型分类及特点Table 4 Classification and characteristics of microscopic constitutive models of hydrate-bearing soils

4 微观力学实验测试技术

当前主要通过三轴剪切及直剪实验对含水合物土力学特性进行研究,尤其在三轴剪切实验方面成果颇丰[125,143-145].微观力学实验方面则可根据其研究方式的不同分为两类,一类是将宏观力学实验与微观测试技术相结合,通过水合物微观结构的演变对宏观力学特性进行分析[106],属于间接方式;另一类是基于三轴剪切等力学实验原理,研发专用于微观力学特性研究的实验装置,直接测量水合物及其界面力学参数的测试技术,如原子力显微镜、微三轴仪等[46,81].

4.1 间接测试技术

间接测试技术一般只测量研究对象的结构,力学性质需要结合三轴实验、直剪实验及声学实验等实验手段.随着现代测试及分析技术的发展,众多具有交叉学科属性的测试技术被广泛应用于水合物研究,而在微观力学领域,则以光学、电学微观仪器为代表.扫描电子显微镜、透射电子显微镜(transmission electron microscope,TEM)、原子力显微镜等能够获取水合物表面微观结构信息,X 射线衍射(Xray diffraction,XRD)、拉曼谱能够实现水合物晶体结构参数分析,计算机断层扫描、核磁共振(nuclear magnetic resonance,NMR) 及成像(magnetic resonance imaging,MRI)则可实现水合物合成及分解过程的结构演化结果,此外也可采用多技术相结合的方式进行综合分析.

4.1.1 计算机断层扫描技术

计算机断层扫描技术被广泛应用于多学科的微观结构测试,其基本原理是:基于结构内不同组分对X 或γ 射线的吸收能力不同,通过入射前后射线强度变化的对比分析,获取结构截面入射方向上的衰减总值,最后利用图像重建技术得到CT 二维图像,借助进一步的衰减系数数字矩阵则可得到三维图像[146].CT 技术的主要优势在于:(1) 技术成熟、分辨率高,能够有效开展孔隙尺度的组分及组构分析;(2) 无损观测,结构完整性高.基于射线扫描,能够有效保证结构对象的完整性;(3) 结合性高,适用性广.能够有效开展基于CT 的联合实验,具有广泛适用性,如CT 三轴剪切实验[147];(4) 分析结果拓展及延续性强.基于CT 扫描的微观孔隙结构表征[148]能够有效用于数值计算等领域研究,如水合物孔隙渗透特性模拟[149].而水合物CT 分析的主要问题为:(1) 水合物与水的区分问题.水合物与水密度相近,因此在CT 图像中两者灰度接近,可通过水中添加溶解盐质[107,150],或使用特殊客体分子[151]实现两相灰度区分.(2) 低含量水合物识别问题.当水合物含量较低时,其对射线吸收较弱,易被其他结构成分掩盖.(3) 水合物动态生长、分解过程中的气水运移及水合物演化识别问题[36].水合物CT 具体应用方面,开展了包括不同水合物孔隙尺度分布[152]、生成及分解过程水合物转化特征[153-154]、渗透特性孔隙结构关系[149]等相关研究,这些成果为水合物宏观特性的微观解释提供了有力支撑.

4.1.2 扫描电子显微镜技术

扫描电子显微镜技术通过聚焦电子束对研究对象进行逐点扫描,成像信号可采用二次电子、背散射电子或吸收电子等,并以二次电子为主.利用SEM 可实时获取研究对象的表面形貌、颗粒关系、孔隙特征等微观信息.在微观组构关系研究方面具有十分重要的作用,尤其适用于颗粒及多孔介质研究对象[155].通过超低温环境实现常压条件下的水合物稳定存在,形成专用于水合物微观结构分析的低温扫描电子显微镜(Cryo-SEM),而依据发射电子的不同可进一步分为钨灯丝型及场发射型扫描电镜.在实际操作过程中,Cryo-SEM 主要难点在于[36]:(1) 由于超低温环境,制样及转移过程中可能发生水蒸气表面凝结,需采取干燥空气等措施避免表面结霜;(2) 喷镀作业时,注意试样表面结构完整性,并保证全程超低温环境;(3) 测试过程避免试样受电子束的损毁破坏,对电子束操作要求较高.SEM 技术在天然气水合物赋存特征[156]、多研究角度下的水合物表面结构特征演化[157]、分解动力学分析[158]等方面得到了广泛应用.

4.1.3 XRD 技术

XRD 技术的原理是基于晶体的周期结构而产生的X 衍射现象,经晶体结构后电磁辐射发生入射方向及强度的变化.被广泛应用于物质结构的定性和定量分析、微观结构分类、完整程度分析等,在化学、力学、医学、生物学、材料学等多学科领域发挥着重要作用.天然气水合物XRD 技术一般通过液氮实现低温,能够有效用于区分水合物晶体类型,是研究水合物晶体结构参数的重要技术手段.水合物XRD 技术实际应用时,要密切关注测试条件度实验结果的影响,尤其是对XRD 图谱的影响.重点关注的测试条件包括:步长选择、扫描速度、累计次数以及测试温度等.目前XRD 技术在水合物晶体结构类型识别[159]、晶体参数分析[160]、生成分解动力学研究[161]以及“自我保护效应”机理分析[162]等方面取得了重要的成果.除上述介绍外,水合物微观研究还包括NMR，MRI，拉曼谱等[163-165]微观测试技术.

微观测试技术与宏观力学实验相结合的宏−微观研究方式有多种,其中CT 三轴应用最为广泛.CT 三轴实验将宏观三轴剪切实验与微观CT 扫描技术相结合,在力学特性分析的基础上结合剪切过程的微观结构观测,实现剪切破坏全过程的宏−微观分析.

4.2 直接测试技术

与间接测试技术需要结合宏观力学实验不同,直接测试技术将微观测试与微观力学实验相统一,能够直接获取微观力学特性及微观结构演变.相对而言,该方式在直接获取微观力学性质方面优势明显,而其主要难点在于力学实验的微观实现.在微观化的过程中,不仅伴随实验装置设计及制造精度的提高,实验操作的难度也将进一步加大.此外一些宏观实验原理在微观层面也将不再适用,这些都极大地限制了直接测试技术的发展,也因此,目前可用于天然气水合物直接测试的微观装置相对较少.

4.2.1 原子力显微镜

原子力显微镜最早由IBM 在20 世纪末所发明,其基础是扫描隧道显微镜(scanning tunneling microscopy,STM),两者则被合称为扫描探针显微镜(scanning probe microscopy,SPM)[166].AFM 将微型测力装置与电子显微镜相结合,能有效开展包括绝缘体在内的各种以固体为研究对象的微观力学实验,在多学科领域得到了广泛应用.天然气水合物研究方面,通过增加温压环境,则可开展各种气体、非气体水合物的微观力学特性研究.

AFM 的结构示意图可见图12,其主要结构包括微悬臂、激光发射器、激光探测器及载物台等[167].其中微悬臂为核心部件,其一端固定另一端则有微小的针尖,即探针.探针的运动由控制器通过调节电流实现控制.激光发射器将发射的激光照射在微悬臂的背面,并反射到由光敏二极管所组成的探测器上,这一发射、接收过程即实现了微小应变的光学信号放大.在样品接触扫描过程中,探针与样品表面相互作用,并引起探针位移(弯曲程度)的变化,进而造成反射光束的偏移,光束偏移蕴含的样品表面信息则被激光探测器所记录.在这一过程中,控制器会自动调节微悬臂的弯曲程度,保持探针与样品表面的排斥力恒定,即两者相对位置始终处于合理范围内.AFM 的工作模式主要有接触模式、非接触模式以及轻敲模式[168](介于接触模式与非接触模式之间),不同模式适用于不同的研究对象,并分别具有不同的优缺点,可根据研究需要进行选取.

图12 AFM 结构示意图[167]Fig.12 Schematic diagram of AFM[167]

对AFM 进行改造,增加温压控制系统,则可形成专用于水合物微观力学特性研究的水合物AFM(如图13 所示),当前,部分学者已经开展了相关研究.例如借助于水合物AFM(二氧化硅微球探针),人们发现微球探针的压入导致了水合物样品相变的发生,并引起显著的塑性变形,探针与水合物之间的黏附力主要由水合物表面似液层及分解液形成的液桥所引起[46].当采用AFM 进行水合物的压痕测量时,应当考虑AFM 尖端引起的热融化效应[169].Huang等[170]则对THF 水合物结晶过程中的微纳米气泡影响过程进行了研究,重点探讨了纳米气泡与水合物间接触角、纳米气泡周围水合物晶体分布模式等演化规律.这些研究结论对于进一步分析水合物微观变形机理、颗粒聚集以及黏附力演化规律等均具有重要的指导意义.

图13 水合物AFM 结构示意图[46]Fig.13 Schematic diagram of hydrate AFM[46]

4.2.2 其他测量仪器

除较为常用的AFM 外,一些具有自主研发性质的“微观机械力学装置”也被用于天然气水合物的微观力学研究.其中,较为广泛的是基于胡克定律而被研发的一种“微机械测力装置”[44],其基本的工作原理如图14 所示.首先在悬臂终端生成水合物颗粒,随后水合物颗粒发生接触与分离,在此过程中记录低弹性常数悬臂终端的颗粒位移,最后根据颗粒位移、悬臂弹性常数,并结合胡克定律则可获得颗粒之间的黏附力.一般此类实验需要进行多组实验,根据实验结果的平均值进行具体分析[171].此外,一些以上述“微机械测力装置”为基础的“改进升级”装置也被相继研发[86],在此不再赘述.

图14 含水合物土的微机械测力装置原理图[44]Fig.14 Schematic diagram of devices to measure force in hydrate-bearing soils[44]

基于宏观三轴实验原理,日本科研人员研发了一种微三轴实验装置[81],可对直径5 mm,高度10 mm 的试样开展三轴剪切实验.同时结合CT 技术,借助于亚毫米及微米尺度图像,实现剪切面、颗粒变形及旋转等局部变形的量化分析.Seol 等[84]也提出了一种全新的CT 结合微三轴实验设备,能够实现成分替换、流体转移、颗粒旋转及破碎的动态观测.

水合物微观力学研究相比宏观性质而言,除了直接的微观力学关系外,还需要对间接的微观细节信息(如组构关系、微观相态分布、动态演化等)进行汇总分析,这些研究内容的有效开展很大程度上依赖于微观力学装置(包括微观测试技术).通过上述总结可知,水合物微观力学研究领域主要集中在近十年范围内,而这一时期也是微观力学装置快速发展的阶段.在宏微观相结合的趋势下,实验装置匮乏与微观力学发展诉求之间的不协调应当引起广大学者们的重视,并作为下一阶段重点解决的问题之一.

5 含水合物土微观力学研究挑战

围绕含水合物土微观力学特性,广大学者开展了从原子(分子)到晶体,从晶体到颗粒(界面),从颗粒到整体(含水合物土)的系统性研究,在取得丰硕研究成果的同时,也面临着诸多挑战.

(1)含水合物土微观结构辨识及定量表征

含水合物土的内部结构性显著,其宏观力学行为响应很大程度上受控于微观结构演化,微观结构的有效测量、准确辨识及定量表征是微观力学特性研究的基础.然而,目前含水合物土微观力学特性研究更多是基于孔隙尺度的内部结构测量、辨识与表征,但是更小尺度的内部结构数据匮乏,导致纯水合物及其与土颗粒界面的力学行为仍不甚清楚.

(2)非均质含水合物土微观力学特性研究

实际自然环境中,水合物生成受多种因素的影响,含水合物土呈现出“非均质性”,决定了其力学特性也是各向异性的.由于制样方法或实验条件的不足,现阶段针对含水合物土的微观力学特性研究,几乎都是以均质水合物为研究对象,如基于均匀制样的力学实验及均匀成核的MD 研究等,迫切需要深入开展非均质含水合物土力学特性研究.

(3)含水合物土多场耦合微观力学研究

实际工程环境中,水合物开采是一个多场耦合的过程,储层变形场与其他物理场之间的耦合作用规律是实现工程地质风险可控的关键.然而,考虑多场耦合作用的含水合物土微观力学实验研究涉及较少,目前仅通过数值模拟的手段实现了力学多场耦合分析,相应的微观力学实验仍然是一个极具挑战性的课题.

(4)微观本构模型构建及拓展研究

本构模型方面,目前已经建立了部分考虑水合物微观孔隙结构以及反映微观力学特性的本构模型,而其适用范围以及参数的确定则是本构研究的主要难点,当前仍缺乏具有广泛适用性同时具有有限模型参数的微观本构模型.另一方面,围绕本构模型的系统性研究不足,基于本构模型的二次开发及其数值模拟研究较为鲜见,这无疑限制了本构模型的进一步应用与推广.

6 结论与建议

天然气水合物的开采过程是多相多物理场的动态演化过程,在相变、渗流、传热等多因素作用下,微观孔隙结构持续变化,最终引起宏观力学特性的演变.因此,天然气水合物宏观力学特性的内在本质是其微观组构关系,微观力学特性研究的重要性不言而喻.整体而言,在宏观力学研究的基础上,水合物的微观力学研究已经进入快车道,同时在宏微观相结合的趋势下,微观、介观到宏观的力学特性统一问题已经成为广大学者研究的焦点之一.而多学科相交叉的综合研究也将进一步推动天然气水合物微观力学研究的快速发展.在本文系统总结天然气水合物微观力学研究现状的基础上,并针对目前研究存在的不足,得到如下几点认识或建议:

(1)微观力学特性测试装置及技术.含水合物土微观力学研究对实验装置及其技术的依赖度较高,目前可用于直接测量微观力学特性的实验装置相对匮乏.此外在多学科交叉背景下,相关学科的最新微观测试技术尤其是光、电、声学测量仪器,应尽快应用到水合物微观力学研究领域.因此,在对现有实验装置进行升级与改造的基础上,应当大力开展水合物微观力学实验装置的研发工作,这也是含水合物土微观力学发展的必由之路.

(2)多种研究方式相结合.含水合物土微观力学研究是系统性体系研究,然而受各种因素及实际操作的限制,不同研究方法之间尚缺乏有效联系的纽带,不同研究方法所得结论相对独立并缺乏对比性.而不同研究方法之间的成果互通无疑是一种有效的研究捷径,同时彼此之间的结果印证也将使得研究结论更具系统性.例如非均匀制样及特殊储层制样是实验研究的一大难题,而根据PFC 研究,通过合理的模型建立及参数选取,能够有效解决制样难题,在定量表征方面也具有一定优势;基于理论分析的本构模型建立可对实验结果进行修正与完善,其二次开发又可实现数值模拟分析.

(3)微观力学特性本构研究.在构建反映微观力学特性本构模型方面,应当避免过多的模型参数.首先对建模框架进行优选,选取相对成熟并且具有广泛适用性的经典本构模型作为建模框架,通过引入微观关联参数构建起反应微观特性的本构模型.而模型参数应尽可能具有实际物理意义,尽量避免“硬性”参数的引入,降低参数确定难度.此外还应开展基于微观本构模型的优选方法研究,尤其是针对不同储层特性、赋存模式以及多因素影响,这需要进一步开展大量的微观力学实验研究.

(4)典型储层及赋存模式研究.与现阶段所取得的微观力学研究成果相比,典型水合物储层研究略显不足.尤其是全球水合物资源的90%以上赋存于泥质粉砂型储层(裂隙型水合物)中,而微观力学研究尚未开展针对泥质粉砂水合物的有效研究,因此应逐步开展泥质粉砂等典型水合物储层的微观力学研究.特别是针对开采过程中的水合物赋存模式动态演化微观力学研究也应当尽快开展.针对上述研究可重点关注DEM 及微观测试技术的应用.