基于二维马尔可夫链的武汉长江公铁隧道地层识别

程利力， 陈 健， 陈 睿， 魏林春

(1. 华中科技大学 土木与水利工程学院， 湖北 武汉 430074； 2. 中国科学院武汉岩土力学研究所， 湖北 武汉 430071； 3. 上海隧道工程有限公司， 上海 200032)

岩土体在长期形成过程中受到地质成分、沉积环境等因素的作用，导致地层分布存在一定的不确定性，即地层变异性[1]；同时地下工程不同于结构工程，这种不确定性很难进行准确分析，传统方法往往根据已有钻孔数据进行简单线性插值从而获取未知区域的地层分布，这种方法简单方便但计算效果有待提高[2]。

目前，一些学者提出了通过概率预测方法对地层变异性进行模拟。谭卓英等[3]采用一种地层识别方法，对钻孔指标进行有效检测。赵磊等[4]基于支持向量机的方法，根据测井数据建立相关地层识别模型。Cao等[5]利用贝叶斯方法对静力触探试验数据进行识别，从而统计了土层的厚度和边界。刘毅等[6]对地质数据采用主成分分析法进行地层识别，从而获取有效的地质信息。Ching等[7]基于静力触探试验数据，采用小波变换识别最大值，从而对土层进行划分。这些研究主要针对单个钻孔进行地层识别，然而对于多钻孔的地层识别研究较少。Qi等[8]提出了一种基于有限钻孔数据的多变量自适应回归样条法，从而可以预测较为合理的地质剖面。Li等[9]基于二维马尔可夫链模型，为研究钻孔之间的地层识别提供了一种较为简便的方法；其中，地层各方向的转移概率是十分重要的模型参数，纵向转移概率矩阵可以由钻孔资料获得，而横向转移概率矩阵受钻孔数据的限制，往往不足以直接估计[10]。

基于上述讨论，为推进地层识别在岩土工程中的应用，本文采用二维马尔可夫链理论模型，间接获取了横向转移概率矩阵，给出了地层识别的算法流程图；最后，以武汉长江公铁隧道武昌岸右线部分钻孔资料为例，进行了有效的地层识别并验证了该方法的可靠性，为以后实际地铁工程施工提供有效的地质参考资料，从而进行相应的工程分析。

1 工程基本概况

1.1 工程概况

武汉长江公铁隧道是国内第一条横贯长江的公铁两用隧道，隧道包括左右两线，单线隧道长2590 m，由武昌工作井始发，经过武昌堤防后穿越长江，到达汉口工作井。图1为武汉长江公铁隧道地理位置图。

图1 武汉长江公铁隧道地理位置

1.2 水文地质条件

武汉长江公铁隧道位于长江河床及两岸一级阶地区。长江江底环境复杂，其中水、泥、沙较多。长江两岸建筑物密集，路网发达，局部有暗塘、暗沟分布，地下管线、排水管道密布，其中盾构段穿越粉质黏土、粉细砂土、泥岩等复合土层，地质环境复杂，属于中等复杂程度，工程重要性程度为一级，因此，地质资料的准确与否对隧道工程施工安全十分重要。图2为武汉长江公铁隧道水文地质图。

图2 武汉长江公铁隧道水文地质

1.3 地质钻孔信息

由于汉口岸一带建筑物密集，地下管线众多，对场地环境扰动过大，因此，本文采用武汉长江公铁隧道武昌岸右线部分地质钻孔资料；按城市轨道交通岩土工程勘察规范[11]，该工程为中等复杂场地，地下工程勘探点间距宜取30～50 m，因此选取区域大小为 186 m×60 m的工程场地；同时考虑到盾构掘进一环约为2 m，故选取大小为2 m×1 m的单元网格进行划分。该工程地质钻孔的相对位置如图3所示，为简便起见，将钻孔编号为Z1～Z6，将各地质钻孔进行投影，得到各钻孔的地层信息剖面如图4所示，该地层主要有6种土体类

图3 钻孔相对位置

图4 各钻孔反映的地层信息

型，其中A为杂填土，B为3-1粉质粘土，C为3-2粉质粘土，D为粉细砂土，E为中粗砂土，F为泥岩。

2 二维马尔可夫链模型

2.1 二维马尔可夫链的基本原理2.1.1 二维马尔可夫链的条件概率

地质沉积过程可以看作马尔可夫过程[12]，利用马尔可夫链的前提条件是土体转移状态符合一阶马尔可夫性；岩土体在横向和纵向上的沉积不同，假设分别由两个不同的一维马尔可夫链组成，然后将其按Elfeki等[13]的方法进行耦合，其在二维空间的分布如图5所示，左右两侧单元可以由钻孔资料获取，按行从左到右依次进行模拟得到未知单元的状态，其中具体步骤见文献[14]。在已知单元(i-1，j)，(i，j-1)，(Nx，j)的状态分别为Zp，Zq，Zw的情况下，单元Zn的条件概率为：

图5 二维马尔可夫链状态序列

p(Si,j=Zn|Si-1,j=Zp,Si,j-1=Sq,SNx,j=Zw)

(1)

2.1.2 纵向转移概率矩阵

纵向转移计数矩阵可由钻孔资料统计获得，记为Cv：

(2)

则纵向转移概率矩阵为Rv为：

(3)

(4)

2.1.3 横向转移概率矩阵

Walther相序定律[15]认为同一剖面地层在各方向的沉积顺序是相同的，但横向的沉积规模一般大于纵向的沉积规模，不能直接采用纵向转移概率代替横向转移概率。通常，地层沉积往往由走向、倾向、倾角三要素表示，而其中地层倾角的作用可量化为横向与纵向延伸长度之比，记为K[16]，由于地层延伸长度在横向远大于纵向，于是在横向上地层向自身发生转移的概率比较大，则横向转移计数矩阵可表示Ch为：

(5)

横向转移概率矩阵Rh为：

(6)

(7)

2.1.4K值估计

地层在横向与纵向延伸长度之比K可以利用地质钻孔进行估计获得，其基本思路类似于极大似然估计，不同K值得到的土体类型不同，将各钻孔区间的转移似然度相乘即为总的土体转移似然度Q；当Q取得最大值时，K值最优。

Q=Q1…QN-1

(8)

Qz=R(Zi,1=Zm,Zj,1=Zo)R(Zi,2=Zm+1,
Zj,2=Zo+1)…R(Zi,Ny=Zn,Zj,Ny=Zp)

(9)

式中：Zm,…,Zn为第z个地质钻孔剖面区间内对应的土体类型；Zo,…,Zp为第z+1个地质钻孔剖面区间内对应的土体类型；Qz为钻孔z与钻孔z-1两侧钻孔之间的土体类型转移发生概率的似然度；Q为钻孔1到钻孔N所有相邻两钻孔之间转移发生概率的似然度相乘之积。图6为各地质钻孔剖面。

图6 地质钻孔剖面

2.2 二维马尔可夫链的计算流程

本文运用MATLAB软件，结合理论公式推导，通过二维马尔可夫链进行模拟，相应的流程如图7所示。

(1)确定建模区域，读取钻孔数据，确定地质剖面土体类型；

(2)利用适当的采样间隔对横向和纵向区域进行离散化，得到马尔可夫链网格；

(3)将地层信息映射到相应的地质单元，计算横向和纵向转移概率矩阵；

(4)估计地质剖面内横向和纵向延伸长度之比K；

(5)对地质剖面进行二维马尔可夫链模拟得到各地质单元的土体类型；

(6)输出模拟得到的地质剖面结果。

3 结果分析与讨论

3.1 转移概率矩阵估计3.1.1 纵向转移概率矩阵

本文针对该工程区域的6个钻孔地质信息，以钻孔Z1，Z3，Z4，Z6作为约束钻孔，钻孔Z2，Z5作为参照钻孔，由式(2)可得纵向转移计数矩阵，如表1所示；将表1的计算结果代入式(3)可得相应的纵向转移概率矩阵，如表2所示。

表1 纵向转移计数矩阵

表2 纵向转移概率矩阵

3.1.2 横向转移概率矩阵

在求解横向转移概率矩阵之前需先确定K值，分别取1，2，…，20等20个不同的K值，由式(5)可得横向转移计数矩阵，如表3所示；其中，K取12时对应的横向转移概率矩阵如表4所示。

表3 横向转移计数矩阵

表4 横向转移概率矩阵(K=12)

利用式(8)(9)依次得到不同K值相应的转移似然度Q，并绘制Q随K的变化曲线，如图8所示。

图8 lg(Q)随K的变化曲线

由图8可知，土体转移概率似然度Q随K值的增大先增大，然后逐渐趋于稳定，其中， 当K为12时，似然度Q最大。

3.2 地层识别结果

采用上述Z1，Z3，Z4，Z6作为约束钻孔进行二维马尔可夫链模拟，以钻孔Z2，Z5作为参照，用于对以下地质剖面实现的准确性评估，运用蒙特卡罗思想[17]，对不同K值进行500次实现；由于K=4与K=20之间似然度比较接近，因此图9给出了K取1，4，8，12，16，20时的一次典型地层模拟结果。

图9 不同K值时的一次典型地层模拟

为描述不同K值局部钻孔Z2，Z5所对应的土体类型效果，采用模拟地层与实际地层分布的误差率来进行定量评价：

(10)

式中：Ni为i号钻孔单元总数；Zi为i号钻孔模拟单元与实际单元的误差个数。

图10显示了预测钻孔不同K值的误差率。由图9，10可知，当K=1，4，8时，在横向上同一土体类型转移概率较小，地层钻孔之间缺乏较好的横向延伸性，因此预测钻孔Z2，Z5的模拟效果不佳，误差率较大；当K=16，20时，地层在横向上同一土体类型发生转移的概率会增大，因此不同土体类型在横向上很难进行转移，这对预测钻孔Z2，Z5的模拟效果也会产生一定影响。由上述结果也揭示了K为12时地层识别结果相对较优。

图10 预测钻孔不同K值的误差率

其中，当K=12时，表5列举了钻孔Z2，Z5处20，40，60 m深度下土体单元分别为不同土体类型的概率值。

由表5可知，该方法将概率最大的土体类型作为该单元的土体类型。同时，根据K=12时的地层模拟结果，表6列举了部分单元位置相应的土体类型，其可作为虚拟钻孔用于工程地质资料分析。

表5 土体单元概率

表6 单元土体类型

3.3 结果验证

利用二维马尔可夫链模拟得到估计钻孔Z2，Z5的地层信息；同时采用钻孔Z1，Z3，Z4，Z6进行简单线性插值得到钻孔Z2，Z5的地层信息，进而对比分析两种方法识别结果的误差率；图11给出了当K=12时，钻孔Z2，Z5模拟结果与实测数据的地层分布情况。

由图11可知，钻孔Z2，Z5利用线性插值法的误差率分别为15%，10%；而利用二维马尔可夫链模拟的误差率为8.3%，6.6%；相比之下，二维马尔可夫链模拟有效降低了地层识别结果的误差率，使得地层识别结果更加接近于真实情况，从而验证了该方法的可靠性。

4 结 论

本文利用二维马尔可夫链理论模型，以武汉长江公铁隧道右线部分钻孔资料为例，进行了有效的地层识别，分析得到以下结论：

(1)在该工程中，基于有限钻孔资料，运用二维马尔可夫链模型进行有效的地层识别，并较简单插值法的误差率有所降低，但各钻孔之间长度不宜过长，否则会对地层识别效果产生影响。

(2)采用二维马尔可夫链模型进行地层识别，充分考虑了地层变异性对地层识别的影响，可为实际工程提供较为准确的地层识别结果，特别是可为弥补地质勘查孔位不足时提供较为详细的工程地质资料。