木里地区永久冻土层的形成及对天然气水合物的成藏影响

王莉平，朱英豪，李 宁，徐拴海，田延哲，刘乃飞，王文丽

(1.西安理工大学 土木建筑工程学院，陕西 西安 710048；2.中煤科工集团西安研究院有限公司，陕西 西安 710077；3.陕西省岩土与地下空间工程重点实验室，陕西 西安 710055)

天然气水合物是由水分子和小分子组成的具有笼状结构的白色或浅灰色冰雪状结晶化合物，因其中气体分子以甲烷(CH4)为主(体积分数＞90%)，也称为甲烷水合物或水合物(Gas Hydrate)[1]。天然气水合物是气体分子和水分子在低温高压条件下形成的白色结晶状物质，主要赋存于海底沉积物、陆域永久冻土带及一些深水湖泊底部沉积物中[2]，由于其能量密度高、分布广、规模大、埋藏浅、成藏物化条件优越，被视为21 世纪最具潜力的代替煤炭、石油和天然气的新型洁净能源之一[3]，2017 年11 月3 日，国务院正式批准将天然气水合物列为新矿种[4]。

我国是世界上第三冻土大国，在青藏高原和大兴安岭地区分布着大片冻土区，多年冻土面积达2.15×106km2，占国土总面积的22.4%[5]。1998 年，据中国科学院政策局有关人员研究，“青藏地区具有形成天然气水合物矿藏的条件和可能，值得引起注意”。1999 年，徐学祖等[6]对青藏高原地区的多年冻土与已发现的赋存有天然气水合物的极地冻土进行对比分析，同时结合几种水合物的标准相图，认为我国青藏高原高山多年冻土有可能赋存有以硫化氢、乙烷和丙烷为主体的重烃类天然气水合物，且其埋藏深度较浅(可能为100～1 000 m)。陈多福[7]、吴青柏[8]等分别定量分析了青藏高原羌塘盆地中可能的水合物储层的顶底界。祝有海等[9]依据青海祁连山木里地区的实测气体组分，结合年平均地温、地温梯度和冻土层厚度等资料，定量分析了天然气水合物形成的热力学条件，认为祁连山木里地区也基本具备天然气水合物的形成条件。2008-2009 年，中国地质调查局组织实施“祁连山冻土区天然气水合物科学钻探工程”，于2008 年11 月5 日在木里地区DK-1 孔井深133.5～135.5 m 处首次发现天然气水合物实物样品，之后分别在11 月7 日和11 月10 日再次发现水合物。2009 年5 月31 日至10 月11 日在DK-2、DK-3 孔中再次钻遇天然气水合物，从而证实我国冻土区存在天然气水合物[10]。

从天然气水合物中烃类气源的角度，一种观点以煤型气为主，认为本区产出的天然气水合物属于偏煤型气的混合气成因[11]。一种观点以油型气为主，认为天然气水合物属于热解成因，主要为石油伴生气或油型气，与煤型气关系不大[12-14]。无论是煤型气，还是油型气，主要的烃源岩分布在上三叠统尕勒得寺组湖相泥岩和侏罗系窑街组。这2 套烃源岩于中侏罗世晚期-早白垩世(165 Ma～100 Ma)进入生烃高峰，是天然气水合物最有利的气源岩，通过青海湖周缘原生砂楔群分析和年代学测试，表明祁连山中更新世早期已进入冰冻圈[15]，即先有气藏，后有永久冻土；此外，还有从气候变化的角度，分析气候变暖永久冻土层变薄之后对天然气水合物赋存的影响[16-17]等。

永久冻土层的形成是青海木里地区天然气水合物成藏的重要条件，目前的相关分析多是基于地质和气候调查取得的间接证据进行判断，缺乏直接定量的分析手段。笔者基于青海木里地区气候、地层的调查事实，分析现有永久冻土层形成时的气候条件及地层温度，利用FLAC3D模拟计算永久冻土层的形成过程、形成时间、稳定底边界及其中的温度梯度，并与调查现状进行对比分析，从定量的角度回答永久冻土层的形成，同时分析此种条件下天然气水合物的赋存条件，并与现状进行对比分析，推测天然气水合物可能的形成过程，相关成果可为木里地区天然气水合物的勘探和开采提供一定的思路。

1 木里地区天然气水合物的勘探现状

近十多年来，我国在青海木里地区天然气水合物的研究方面取得了巨大的进展，研究主要集中在钻探调研、成分分析、来源推断及气候变化可能对天然气水合物赋存带来的影响方面。如中国地质调查局先后在木里地区钻探DK-1-DK-12 共12 口调查井，其中DK-8 兼为试采井，在DK-1、DK-2、DK-3、DK-7、DK-8、DK-9、DK-12 中发现赋存的天然气水合物，在DK-4、DK-5、DK-6、DK-11 中发现与天然气水合物相关的异常现象，在DK-10 中发现异常高压浅层天然气层[10,14,18-19]。神华青海能源开发有限责任公司项目“青海省天峻县聚乎更煤矿区三露天天然气水合物调查评价”共组织施工13 口天然气水合物探查钻孔，其中4 口井探获天然气水合物实物样品，分别为DK8-19、DK11-14、DK12-13 和DK13-11[20]。以上各钻孔的分布如图1 所示。

图1 木里煤矿天然气水合物研究区地质及钻井分布示意图(据参考文献[21]修改)Fig.1 Geological condition and boreholes distribution of the natural gas hydrate study area in Muli Mine(According to Reference [21],modified)

2 木里地区多年冻土的现状及形成分析

2.1 基本特征

青海木里地区位于祁连山中部，是介于托来山和大通山之间的一个山间盆地，高程一般为4 000～4 300 m，是祁连山冻土区的核心，除局部融区外，多年冻土连续分布，年平均地表地温最低-2.4℃，实测冻土层厚度60～95 m，并常见厚层地下冰。平均地表地温和多年冻土厚度实测数据[5]见表1。

表1 木里地区多年冻土区平均地表地温和厚度实测数据[5]Table 1 Average ground surface temperatures and thicknesses of permafrost in Muli Region[5]

表1 中的数据表明，木里地区多年冻土区的年平均地表地温为-1.95℃，多年冻土平均厚度为71 m，假设地温沿着冻土层深度方向线性变化，则地温梯度为0.027 3℃/m，即2.73℃/hm。

王超群等[22]对木里地区DK-9 孔进行了连续3 a的地温监测，监测结果表明多年冻土厚度约为160 m，年平均地温为-3.3℃，多年冻土上限为-2.5 m。冻土层内的地温梯度为1.38℃/hm，冻土层以下的地温梯度为4.85℃/hm。金春爽等[23]在DK-1 钻孔进行天然气水合物勘探时，钻探过程中同时开展井内温度测量，提钻后24 h(使钻探的热效应对岩层影响最小)分别对3个充满泥浆的钻孔进行温度测井。从钻孔泥浆温度测量结果来看，DK-1 孔中，钻孔泥浆上部约20 m 温度不断降低，可能反映了大气与冻土层相互作用的结果；在20～115 m，泥浆温度不断增大，呈现出一定地温梯度的变化特征，其变化斜率代表着冻土层内地温梯度值，约为1.39℃/hm(R2=0.95)；在120 m 深度以下，泥浆温度增大的斜率与115 m 深度以上明显不同，为2.86℃/hm(R2=0.97)，冻土底界出现在115～120 m 深度处。

木里煤田聚乎更矿区呈NWW-SEE 向展布，东西长约19 km，南北平均宽约4 km，面积约76 km2，总体上为一复式背斜构造，由1 个大背斜和2 个小向斜组成。其中北向斜分布有三井田、二井田和一露天3 个井田，南向斜由四井田、一井田、三露天和二露天组成，如图2 所示。

图2 木里煤田聚乎更矿区井田分布[24]Fig.2 Wellfield distribution in Jvhugeng Mine,Muli Coalfield[24]

依据中铁资源集团海西煤业聚乎更矿区四井田露天矿(首采区)《边坡工程地质勘察与稳定性评价报告》(2013 年)[25]的现场勘察，四井田矿区内多年冻土分布很不均匀，东部高程较低，多年冻土的厚度较小，为50～70 m，西部高程较高，多年冻土的下限最深超过150 m。在四井田北部非工作帮设bk2-4 孔为地温监测孔，监测结果表明该区多年冻土上限为-5.75 m，多年冻土厚度为54 m 左右，多年冻土中的温度梯度约为2.85℃/hm，多年冻土之下的温度梯度约为4.36℃/hm。

由以上数据可以看出，木里地区多年冻土层的厚度差异很大，在50～160 m 之间变化，相较而言，发现天然气水合物的三露天井田中多年冻土的厚度更大，DK-1、DK-9 孔揭示的多年冻土层底界分别出现在115～120 m 和160 m 深度处。另外，多年冻土层内的温度梯度在1.3～3.0℃/hm，永久冻土层越厚，其中的地温梯度越小。

2.2 木里地区多年冻土层形成时间分析

胡道功等[15]通过青海湖周缘原生砂楔群分析和年代学测试，揭示祁连山现代冻土为形成于中更新世早期倒数第3 次冰期(770 ka)、倒数第2 次冰期(40 ka)和末次冰期(13 ka)的冻土演化而来，确定各冰期时期年均气温为-10～-5℃。金会军等[26]对2 万年来的中国多年冻土形成演化进行了分析，大致分为7 个阶段：晚更新世LGM(20 000～10 800 a BP)多年冻土强烈扩展，达到LPMax；早全新世气候剧变期(10 800 至8 500～7 000 a BP)多年冻土较稳定但相对缩减阶段；中全新世HMP(8 500～7 000 至4 000～3 000 a BP)多年冻土强烈退化阶段，多年冻土缩减到LPMin；晚全新世新冰期(4 000～3 000 至1 000 a BP)冻土扩展阶段；晚全新世中世纪暖期(1 000～500 a)多年冻土相对退化阶段；晚全新世小冰期(LIA，500～100 a BP)冻土相对扩展阶段，以及近代升温期(近百年来)多年冻土持续退化阶段。

笔者团队于2013 年7 月在木里地区聚乎更矿区四井田露天矿采坑一标段的南部工作帮中段采取冰层试样(见冰深度15～18 m)进行了14C 检测，测龄分别为2 136 a 和1 348 a。

结合以上分析及现场实测结果，初步有以下结论：

(1) 木里地区自进入第四纪以来，多年冻土层历经变迁：生长-融化-生长-融化-···，并非在固定时期形成并保持至今；

(2) 在中全新世HMP(8 500～7 000 至4 000～3 000 a BP)，木里地区的多年冻土可能融化殆尽；

(3) 木里地区目前的多年冻土可能形成于晚全新世新冰期(4 000～3 000 至1 000 a BP)的冻土扩展阶段。

因此，在之后的多年冻土形成的模拟计算中，拟以晚全新世新冰期(4 000～3 000 至1 000 a BP)的气候作为温度边界条件。

2.3 多年冻土层形成时的气候条件

全新世晚期(4 000～3 000 a BP)气候开始显著变冷。祁连山的敦德冰心记录在4 000 a BP 开始降温，在2 800～2 700 a BP 达到极低值。自此，温度波动下降直至1 000 a BP 结束，此时段是全新世较寒冷的新冰期阶段。通过对比青藏高原的多年古冻土和古冰缘现象，推算当时高原多年冻土下界比现今普遍低约300 m，当时年平均气温比现在低约2℃[26]。

3 木里地区典型烃类气体的相变曲线

截至目前，祝有海等[10]在DK-1 井100 m 深处经排水集气法收集气体及泥浆气共3 个；刘昌岭等[27]在DK-2 井140～480 m 内，采用真空顶空法收集样品中的水合物分解气6 个，在DK-3 井142 m 和395 m，采用真空顶空法收集样品中的水合物分解气2 个；黄霞等[14]采集了DK-2、DK-5、DK-6 井中的岩心游离气共59 罐。相较而言，刘昌岭等[27]采集的DK-2 井中的样品为水合物分解气，比岩心游离气更能反映天然气水合物的构成，且DK-2 井中的数据更为翔实，因此，本文以DK-2 孔中的天然气水合物储层状况作为模拟及分析对象，DK-2 井中天然气水合物的气体组成比例见表2。

表2 DK-2 钻井中烃类气体组分摩尔比[26]Table 2 Molar ratio of hydrocarbon gas components in DK-2 borehole[26]

黄霞等[13]对木里地区采集到的天然气水合物进行了拉曼光谱检测，结果表明木里地区的天然气水合物属于II 型水合物，即除甲烷外，还含有较高的乙烷、丙烷等重烃成分。表2 中的烃类气体样品都归属为II 型水合物，但可以看出样品1 与样品2-样品6 的气体组成差异显著：样品1 中的甲烷摩尔比为34.85%，样品2-样品6 中的在60%左右，样品1 中重烃含量远高于样品2-样品6 的。这些组分上的差别相应地也导致了成藏时温压条件的差异，为了更有针对性地进行水合物成藏研究，特地根据气体组分摩尔比的不同，将天然气水合物划分为A 类和B 类。可以看出，对DK-2 井而言，A 类天然气水合物埋藏深度约在井深149 m，B 类天然气水合物的埋藏深度在井深253 m 以下，取样品2-样品6 的平均值代表B 类水合物。

利用Sloan 的CSMHYD 软件对A 类和B 类天然气水合物形成的温压条件进行计算，得到的结果如图3 所示。

由图3 可以看出，当温度小于2℃时，A 类和B 类天然气水合物成藏所需的压力基本相同；当温度大于2℃时，A 类天然气水合物成藏所需的压力逐渐超过B 类天然气水合物，且二者之间的差值越来越大；当存在永久冻土层时，即地层温度小于0℃时，其中A 类和B 类天然气水合物成藏所需的压力在250～500 kPa，即便只考虑土体自重，此压力也非常易于获得，在现有地层情况下，埋深大于13～25 m 即可达到成藏 条件。

图3 A 类和B 类天然气水合物形成时的温度压力条件Fig.3 Temperature and pressure conditions for the formation of Class A and Class B gas hydrates

4 永久冻土层形成过程计算

4.1 地层建模及相关参数

利用有限差分程序FLAC3D进行永久冻土层的形成计算。考虑到木里地区钻孔揭示的最深永久冻土层底界为160 m，本次计算地层深度(Z方向)取200 m，宽度(X方向)取10 m，厚度(Y方向)取1 m。

DK-2 孔具体的岩性厚度比例组成[28]见表3 。

表3 DK-2 钻孔中不同岩性的厚度比例Table 3 Thickness ratio of different lithologies in DK-2 borehole

根据以上钻探获得的DK-2 孔岩性分布情况，结合钻孔岩心柱状图的0～200 m 岩层条件，在进行地层建模时，将DK-2 孔0～200 m 地层视为基岩(泥质粉砂岩+粉砂岩+粗砂岩)、水及冰的混合体，混合体中各成分的基本物理参数见表4，其中基岩的指标采用泥质粉砂岩、粉砂岩和粗砂岩[29]的平均值。

表4 地层混合体中各成分物理参数Table 4 List of physical parameters of each component in the stratigraphic mixture

三露天井田天然气水合物钻孔岩心裂隙十分发育，裂隙形式复杂多样[20]，因此，在考虑基岩+裂隙混合体的热传导系数和体积比热容时，依据混合物理论[30]，采用指数加权模型描述混合体等效热传导系数λe，其随冻结温度的变化规律为：

式中：λs、λl、λi分别为岩石基质、水以及冰的热传导系数；wu为裂隙中未冻水含量；n为混合体中的裂隙率。

对于混合体而言，已冻区和未冻区的热平衡方程的差别主要体现在相变项上，相变潜热在等效比热容中予以考虑。混合体的等效体积比热容cv可表示为：

式中：ρs、ρi、ρw分别为基岩、冰及水的密度；cs、ci、cw分别为基岩、冰及水的体积比热容；l为水冰相变潜热；θ为水的温度。

本文将裂隙中的水冰相变限制在-3～0℃内，冻结终了时有wu=0，则裂隙水的未冻水含量方程可表示为：

4.2 温度边界条件

DK-2 孔缺乏直接的温度监测资料，依据附近DK-9 孔长期的地温监测[22]：地表处平均温度为-3.3℃，以正弦曲线模拟1 a 内的温度变化，振幅采用月平均最高温/最低温度与年平均温度的差值，为11.45℃。

由2.2 节的分析可知，晚全新世新冰期(4 000～3 000 至1 000 a BP)可能是目前木里地区多年冻土的形成及成长阶段，此阶段的气温较现今低约2℃，因此，将以上的正弦曲线下移2℃，即年平均地表温度为-5.3℃，振幅保持不变，以此作为气温边界条件。

式中：以某年1 月1 日零点为起算点，s为距离起算点的时长，s；θ1为s时间对应的温度，℃；L为一年的总时间，则L=365×24×3 600，s；b1、c1为地表的平均温度和温度振幅，分别取-5.3℃和11.45℃。

模型两侧为绝热边界，模型底部边界取热流边界，热流密度为0.02 W/m2。

4.3 地层初始温度场

依据2.2 节的分析，木里地区的多年冻土在中全新世HMP(8 500～7 000 至4 000～3 000 a BP)时可能融化殆尽，成为季节性冻土。根据施雅风[31]成果，中国HMP 出现于8 500～3 000 a BP，其中稳定暖湿的鼎盛阶段在7 200～6 000 a BP，青藏高原南部气温比现今高4～5℃，则依据DK-9 孔的地温监测，取此时的地表处平均温度为-3.3℃+5℃=1.7℃，幅值仍取11.45℃。不同埋深处的初始地温按如下公式进行估算[32]：

式中：z为距离地表的深度，取绝对值，m；t为距离起算日期的时间，d；p为振动周期，可取365 d。

地温振幅≤ ±0.1℃的深度称作地温年变化深度，该处的地温称作年平均地温，年变化深度以下的区段内，地温基本上不随时间变化，可视为稳定温度场[33]。由式(7)计算可得：地温年变化深度为16 m，年平均地温为1.7℃，该深度以下的地温场梯度按照3℃/hm考虑。

5 计算结果及分析

5.1 永久冻土层形成时间及底边界深度

随着计算的进行，永久冻土层的底边界不断下移，当计算时间持续到170 a 及以后时，永久冻土的底部边界基本保持不变，深度维持在130 m 左右，具体计算结果如图4 所示。

图4 永久冻土层形成过程中地层中的温度分布云图Fig.4 Cloud map of the temperature distribution in the stratum during permafrost formation

由2.1 节分析可知，DK-2 钻孔缺乏直接的冻土层勘察资料，其附近钻孔DK-1 的永久冻土底边界深度为115～120 m，DK-9 钻孔揭示的冻土底边界为160 m左右，这2 个钻孔，尤其是DK-1 钻孔，与本文的模拟对象DK-2 钻孔位置很近，间接证实了永久冻土模拟计算具备一定的可靠性。

5.2 永久冻土层中的温度及温度梯度

从25 m 左右的深度开始，地温曲线呈现线性变化规律(θ=0.016 4d-2.199 4，R2=0.996 7)，温度梯度为1.64℃/hm，符合2.1 节中分析的此区域内永久冻土层中的温度梯度，且由于永久冻土层较厚，温度梯度偏小，具体计算结果如图5 所示。

图5 永久冻土层底边界稳定时地层温度与深度关系Fig.5 Relationship between formation temperature and depth when the bottom boundary of the permafrost layer is stable

5.3 天然气水合物形成时的温压条件

假设地层中仅有自重压力，将模拟计算得出的地层自重压力及对应温度与A 类和B 类天然气水合物形成时的温压曲线绘制在同一张图中，如图6 所示。

由图6 可以看出，当地层埋深超过20 m 时，即地层压力大于300 kPa 时，永久冻土层中即具备了形成A 类和B 类天然气化合物的条件。目前钻孔揭示出来的A 类和B 类的储层分别在149 m 左右和大于250 m 处，可能的原因有2 种，一种为深部地层中烃类气体迁移时仅到达140 m 左右处，可能由于地质构造方面的原因导致气体无法进一步向上迁移，无论地层处于冻结状态还是融化状态，天然气水合物仅随之改变存在状态(固态或气态)。第二种原因是，烃类气体的确迁移至地层埋深较浅的位置，但由于在第四纪时期，气候历经变迁，永久冻土层并非处于持续稳定的状态，而是融化-生长状态交替出现，使得永久冻土时期形成的天然气水合物在气候转暖时又动态分解为烃类气体，并沿着裂隙向上迁移释放至大气中去。

图6 永久冻土层底边界稳定时地层内温压曲线与A/B 类天然气水合物形成时的温压曲线Fig.6 Temperature-pressure curves in the formation when the bottom boundary of the permafrost is stable and the temperature-pressure curves during the formation of class A/B natural gas hydrate

6 结 论

a.青海木里地区天然气水合物成藏区域现存永久冻土可能形成于晚全新世新冰期(4 000～3 000 至1 000 a BP)的冻土扩展阶段，计算结果表明，当晚全新世新冰期降温时，历经约170 a，永久冻土的底边界基本达到稳定，深度约为130 m，永久冻土层中的温度梯度约为1.64℃/hm，与现状勘察的结果相当吻合。

b.当永久冻土层形成时，在深度20 m 及以下即具备A 类与B 类天然气水合物的形成条件，A、B 类天然气水合物均在140 m 深度以下，可能的原因是：地质构造使得深部烃类气体只能迁移至140 m 深度以下，随着永久冻土层的形成，烃类气体在适当的温压作用下形成天然气水合物；二是烃类气体可能迁移至更浅层的地层中，但由于多年冻土的反复演化，使得永久冻土时期形成的天然气水合物在气候转暖时发生动态分解，以气态形式进入大气。

c.与北极地区相比，木里地区天然气水合物分布相对分散，可能的原因有气候作用下多年冻土的反复演化，使冻土处在冻结-融化-冻结的交替变化中，由此导致地层中的温度和压力始终处于变动中，进一步导致天然气水合物储存的分散。

d.建议在未来的天然气水合物探矿勘察中重点关注永久冻土层变化不大的地区，深部地层同时具备烃类气体的形成条件，如羌塘盆地等。