2020, Vol. 52
2. 长崎大学 工学院，日本 长崎 8528521;
3. 北京交通大学 土木建筑工程学院，北京 100044;
4. 城市地下工程教育部重点实验室(北京交通大学), 北京 100044
2. Graduate School of Engineering, Nagasaki University, Nagasaki 8528521, Japan;
3. School of Civil Engineering, Beijing Jiaotong University, Beijing 100044, China;
4. Key Laboratory of Urban Underground Engineering (Beijing Jiaotong University), Ministry of Education, Beijing 100044, China
一般认为，与地面建筑相比，隧道等地下结构在地震时由于受到围岩约束，震动幅度小且震动次数少，抗震性较好.因此，以往关于隧道的抗震设计相对较少.然而，近几十年来发生的隧道震害表明了该认识具有片面性，例如，1995年日本阪神地震、2008年中国四川汶川地震以及2016年日本熊本地震等.震区内大量隧道受到了不同程度的破坏：阪神地震的震害隧道主要集中在20世纪60年代建的浅埋隧道，比如山阳新干线的六甲隧道，其拱顶产生了横向剪切裂缝，另外拱墙交接处有压缩性裂缝[1]；2008年汶川地震震区中隧道震害，除衬砌裂缝外，比较突出的还有隧道洞口部的坍塌[2]；而2016年日本熊本地震，隧道的破坏形式又有新的特点，环形裂缝均匀分布[3].以上实例表明，地震导致的隧道破坏受多种因素影响，机理相对复杂.研究隧道等地下结构受地震力作用下的破坏特征和规律，对其防震减灾有重要的现实意义.
从细部研究隧道的震害，着重于力学方面的探讨，是隧道震害研究的一个重要方面.但因地质条件的复杂性以及地震波的不确定性，这种方法受到很大制约.受地震作用影响下，地表变形运动是一种直接有效的地震作用机制观测方法.有研究表明，2016年熊本地震中，地面震动和变形可能会对长大型结构(如岩体)造成巨大的影响[4]，然而该研究没有给出明确的研究结果.
为此，拟从地震前后地表变形角度，探讨分析隧道的破坏特征，探求其两者之间的相关信息，以便通过地表变形特征划分地下结构的震害影响区域，用于隧道等地下工程的防灾减灾.
目前，分析地表变形时常用的数据是高清晰度垂直和倾斜航空照片以及机载激光测量数据(LiDAR)[5-8].机载激光测量系统由集成航空定位系统和机载激光扫描器构成，其测量数据可以过滤植被影响，最近在岩土工程中的应用越来越广泛[9-13].目前为止，山岭隧道的地震破坏研究中，相应的定量分析研究还很少，且地表运动与地下结构变形之间的相关性鲜有研究.本文通过实例分析，考察隧道周边地表变形运动分布特征，对比其与俵山隧道震害的空间分布特征，探讨地表变形与隧道震害之间的相关性.
1 2016年熊本地震简介2016年4月14日21点26分(日本标准时间，东九区时间)，日本熊本县熊本地区发生震级6.5(Mj)的逆断层型地震，震源深度约为11 km (2016年熊本地震前震)，地震源于日奈久断层的错位.2016年4月16日凌晨1点25分，由于布田川断层发生走滑型破裂，该地区再次发生震级7.3(Mj)地震，震源深度约为12 km(2016年熊本地震主震).日本气象厅将该一系列地震命名为“平成28年(2016年)熊本地震”.图 1为2016年熊本地震各地震震源分布图(数据来源：日本国土交通省气象厅，2016).
|
图 1 2016年熊本地震各地震震源分布 Fig. 1 Epicenter distribution map of the 2016 Kumamoto earthquake |
图 2为主震发生时熊本西区春日观测点加速度波形，该地震波形分为南北方向(NS)、东西方向(EW)和上下方向(UD)3个分量.由图 2可知，该观测点地震加速度波形水平分量(NS和EW分量)大于上下分量，NS分量峰值为605.9 cm/s2，EW分量峰值为551.5 cm/s2，而上下分量峰值仅为405.2 cm/s2；另外该监测点NS分量大于EW分量，表明该区域地表运动方向以南北向运动为主.表 1给出了熊本震区内多个观测点最大加速度值的水平和垂直分量，可以看出，水平分量大于上下分量(南阿苏中松观测点除外)，水平分量中南北分量和东西分量的比重所占不同，表明地震波由于地形条件差异等因素的影响，在各个区域的传播方向会有所不同，如春日观测点区域以南北向传播为主，而上天草市大矢野町区域以东西向传播为主.
|
图 2 主震发生时熊本西区春日观测点加速度波形 Fig. 2 Acceleration of mainshock wave in KASUGA when the mainshock of Kumamoto earthquake occurred |
|
表 1 熊本地区各观测点最大加速度值(2016年4月16日主震) Tab. 1 Measured maximum ground acceleration at different observation sites induced by themainshock of Kumamoto earthquake (April 16, 2016) |
本次前震与主震最大烈度皆为烈度7，是日本自1949年设立烈度7以来，继1995年阪神大地震、2004年新泻县中越地震和2011年311大地震后第4，5次观测到烈度7的地震，也是日本自1949年来第1次数日内在同一区域连续两次观测到烈度7的地震.
2 俵山隧道概况及震害调查俵山隧道位于熊本高森线(县道)上，全长2.057 km，呈东西走向，是连接西原村和南阿苏村的交通要道.该隧道距主震震中直线距离为22.4 km，如图 3所示.在2016年熊本地震中遭受严重损坏.
|
图 3 俵山隧道概况 Fig. 3 Situation of Tawarayama tunnel |
地震发生后，在日本国土交通省九州分局熊本水川国道事务所的帮助下，笔者和研究小组率先进入震区，对隧道震害进行了详细调查，取得了第一手现场资料.调查发现，该隧道的基本地震破坏特征可分为5种模式：衬砌裂缝(图 4(a))，施工缝损伤(图 4(b))，地下水渗漏(图 4(c))，混凝土衬砌剥落、坍塌(图 4(d))以及路面损坏(图 4(b)).在详细调查之后，绘制了隧道损坏的全景视图以进行进一步评估.图 5为S157-S190跨的隧道震害破损全景图(此处S指代隧道跨度).
|
图 4 俵山隧道震害类型 Fig. 4 Classification of the seismic damages to Tawarayama tunnel |
|
图 5 俵山隧道的震害破损全景图(S157-S190)[3] Fig. 5 Panorama view of the seismic damages to the Tawarayama tunnel (S157-S190)[3] |
熊本地震发生后，日本政府和很多民间机构对灾区进行了基础测绘，通过高精度的航空拍摄和卫星遥感等技术，取得了地震后的基本地形数据.由于机载激光测量技术的高穿透性，该技术能够过滤地表覆盖植被的影响.因此，本研究中采用机载激光测量数据.
分析隧道周边地表变形时，利用的数据是2016年地震发生前和地震发生后测量的数字高程模型(digital elevation model，简称DEM).震前的数据由日本国土地理院提供，测量时间为2013年1月—2月，栅格为1.0 m；震后的数据来自日本林业局，测量的时间为2016年4月—7月，栅格为0.5 m.为方便论述，这里把2013年的数据称为震前DEM，2016年的称为震后DEM.利用栅格为1.0 m的DEM，通过合适的计算方法，地表变形计算精度可达0.1 m.
采用组合分类迭代最近点(CCICP)算法计算地震前后的地形变动特征.该算法最初是为了标定移动测绘系统(mobile mapping system，MMS)的点阵而开发[15].其基本原理为相同区域通过MMS扫描多次时，捕捉有差异点集合，分析处理点阵的变化，可得到地形的变化.文中将震前DEM设为S，震后DEM设为T，地表变形的计算过程如下：
1) 将S和T按照等间隔间距划分为小的点集合数据区域.
2) 根据Demantke等提出的主成分分析法(principal component analysis，简称PCA法)[16-17]，对每个点集合进行分类.其基本原理为：首先根据点所处位置以及周围点的分布，计算出各个点的特征值λ1, λ2和λ3，按照特征将点分为线性点集，平面点集和散点集；然后通过特征值的协方差矩阵对点的位置进行评价，如果点集之间是线-面关系或散点-面关系时，则认为二者坐标系不匹配，如果点集之间是点点关系或者面-面关系，则认为二者坐标系匹配，可进行下一步计算.
3) 计算地震前后点-点集合和面-面集合之间的最小距离.点-点距离DPT_PT以及点-面距离DPT_PL的计算公式如下[15-17]：
| $ {D_{{\rm{PT}}\_{\rm{PT}}}} = \sqrt {{{\left| {\left( {\mathit{\boldsymbol{T}} \cdot {p_{\rm{s}}} - {p_{\rm{t}}}} \right)} \right|}^2}} , $ | (1) |
| $ {D_{{\rm{PT\_PL}}}} = \sqrt {{{\left| {\left( {\mathit{\boldsymbol{T}} \cdot {p_{\rm{s}}} - {p_{\rm{t}}}} \right) \cdot {\mathit{\boldsymbol{n}}_{\rm{t}}}} \right|}^2}} , $ | (2) |
式中：T为齐次坐标系中的变换矩阵；ps为震前点集合，pt为震后点集合；nt为PCA方法计算后pt的法向向量.T的表达式为
| $ \mathit{\boldsymbol{T}} = \left( {\begin{array}{*{20}{c}} {{\alpha _{11}}}&{{\alpha _{12}}}&{{\alpha _{13}}}&{{t_x}}\\ {{\alpha _{21}}}&{{\alpha _{22}}}&{{\alpha _{23}}}&{{t_y}}\\ {{\alpha _{31}}}&{{\alpha _{32}}}&{{\alpha _{33}}}&{{t_z}}\\ 0&0&0&1 \end{array}} \right). $ | (3) |
式中：αij(i, j=1, 2, 3)为旋转量，tx, ty, tz为平移量.
4) 通过反复迭代S和T中对应点坐标系之间的距离平方差，获取S和T之间平方差之和的最小值，以此作为S和T的配准，最终获取地震前后地表变形三维分布特征.
4 结果和讨论 4.1 地震前后的地形变化图 6为通过CCICP算法计算得到的地表水平变形位移.可以看出，俵山隧道东侧地表变形主要是以东北-西南向为主，而隧道西侧地表变形主要是以东-西向为主.同时，将本文计算结果与前人计算结果作对比以验证该研究方法的有效性.图 7为Mukoyama等采用3D-GIV(geomorphic image velocity)方法计算得到的隧道周边位移分布[18].比较二者之间的分布特征可以发现两种研究方法的结果基本吻合.
|
蓝色为CCICP计算水平位移值，红色为日本国土地理院(GEONET)公布的电子基准值，单位：m 图 6 地表水平变形分布CCICP计算值 Fig. 6 Horizontal displacement calculated by CCICP algorithm |
|
黑色箭头代表水平位移，Mukoyama et al., 2017年[18] 图 7 隧道周边的位移矢量 Fig. 7 Displacement vector around tunnel |
此外，为了进一步验证该研究方法的准确性，将附近的长阳(Chouyou)基准观测点(编号为EL04930274901，经纬度坐标为32.870N, 130.9962E)作为研究对象.图 6同时给出了地震前后长阳(Chouyou)基准点附近的水平位移分布.蓝色为本研究的计算值，红色为日本国土地理院公布的观测值.地震后，基准点向西南方向移动了0.97 m.该基准点周边6个计算点水平位移的平均值也为0.97 m，由此可以看出本研究方法计算得到的地表位移值和观测值基本一致，从而验证了该研究方法的有效性和准确性.
4.2 地形变形特征与隧道震害之间相关性分析图 8为俵山隧道震害分布与地表垂直变形对应关系.正值代表地表隆起变形，负值代表地表沉降变形.为便于观测，沿着东西方向设置等间隔监测线，由西向东命名为L_1，L_2, …, L_21.由图 8(b)可知，地震导致的地表隆起和沉降比较小.在隧道的西侧入口，地表的沉降为0.57 m，隧道上方其他部位的隆起和沉降比较平缓，最大值为0.5 m.由表 1可以看出，各个观测点加速度的垂直分量远远小于水平分量，由此也可以说明地表的垂直变形较小.
|
图 8 俵山隧道震害分布与地表垂直变形对应关系 Fig. 8 Relation between ground vertical displacement and seismic damages of the Tawarayama Tunnel |
隧道震害调查结果表明，隧道衬砌结构由于竖直方向作用导致的震害较少(图 7(a)).由此可以看出，该案例中地表垂直变形与隧道衬砌结构震害之间的相关性较弱.
图 9为俵山隧道震害分布与地表水平变形对应关系.由图 9(b)可以看出，位于隧道西侧入口(L_1~L_6)上方地表变形的水平位移较大，平均值为1.15 m；位于隧道中部和东部(L_7~L_21)上方的地形，水平位移值逐渐减小，平均值为0.88 m.现场调查结果表明，隧道西侧(S001~S050，对应于地表检测线L1~L6)，由于水平挤压、张拉或剪切作用，隧道衬砌开裂、剥落现象严重，隧道衬砌施工缝、路面也出现大量破损.由此猜想地表水平变形特征能够在一定程度上反映地下结构变形特征和震害分布特征，两者之间的相关信息便于划分地下结构的震害影响区域，用于隧道的防灾减灾.
|
图 9 俵山隧道震害分布与地表水平变形对应关系 Fig. 9 Relation between ground horizontal displacement and seismic damages of the Tawarayama Tunnel |
为探讨地表变形和隧道震害之间的相关性，沿俵山隧道轴线方向作两条辅助线：北侧测线(North_detection_line，简称N_L)和南端线(South_detection_line，简称S_L)，如图 9(b)所示.N_L测线上监测点采用红色方块标记，S_L测线上监测点采用黑色实心圆标记.N_L测线和S_L测线上各个监测点的南北向和东西向位移分量如图 10(a)，(b)所示.此外，定义地表水平位移与隧道轴线夹角的锐角(< 90°)为θ.图 10(c)为两条测线上监测点水平位移矢量与隧道轴向夹角θ沿隧道轴向的分布情况.
由图 10可知，在隧道西端的地表变形中(L_1~L_6)，东西方向的水平位移大于南北方向的值，地表水平运动方向与隧道轴向夹角小于10°.隧道震害调查结果表明，该区域隧道衬砌呈现较明显的东西向轴向拉压变形破坏特征，如图 4(b)所示的跨度S012与跨度S013之间衬砌施工缝压缩破坏和路面受压隆起破坏，图 11所示的跨度S001与跨度S002之间衬砌施工缝张拉破坏.同时，隧道洞口段坑道向西移位10 cm，并且由于挤压作用导致隧道西侧出口处路面隆起15 cm，如图 12所示.除此之外，如图 9(b)所示，在L_1和L_2处，发现了往北0.19和0.13 m的位移，而其他部分的位移则都往南，地形南北向运动方向发生改变，容易导致地层发生错位.隧道震害调查结果发现，在该地形位移南北方向转变位置的下方，隧道衬砌在施工缝处发生南北向错位，如图 13所示.可以看出，强烈的地表水平变形运动特征能够在一定程度上反映地下结构震害的真实情况，如衬砌、路面和施工缝等各种结构的破坏特征，这将有助于划分地下结构的震害影响区域，用于隧道的防灾减灾.
|
图 10 北侧测线(N_L)和南侧测线(S_L)上的位移分量和隧道轴线夹角θ的对应关系 Fig. 10 Displacement components on both detection lines and variation of angle θ along the axis of the tunnel |
|
图 11 施工缝拉伸破坏 Fig. 11 Construction joint damage caused by tension |
|
图 12 隧道出口处的路面破损(西原村侧) Fig. 12 Portal pavement damage of the Tawarayama tunnel (Nishihara Village side) |
|
图 13 跨度S001和跨度S002之间的施工缝处错动 Fig. 13 Dislocation of the construction joint between S001 and S002 |
由图 9可以看出，隧道中东部上方的地形(L_7~L_21, 南阿苏村方向)，东西向水平位移逐渐减小，南北向水平位移逐渐增大，且位移值超过东西向位移值，θ表现出增大趋势，表明该区域地表变形方向为西南-东北方向.表 1中南阿苏村中松观测点位于隧道东部区域的东北侧，该观测点的南北向峰值加速度值大于东西向峰值加速度值，这与地表变形特征吻合，从而表明地表变形运动的方向能够在一定条件下(如排除地震断层等特定地质构造的影响等条件)反映地震波的传播方向.
地表西南-东北向运动特征表明地震波的传播方向与隧道轴线方向斜交(如图 10(c)所示，θ介于20°和90°)，与隧道东西向轴线方向斜交.该结论验证了文献[3]关于地震波传播方向的猜想.现场调查结果表明，虽然该区域地表的水平变形远远小于隧道西端，但由于衬砌质量、岩体的地质条件等因素，该区域隧道破坏仍然较为严重.如沿着隧道轴向约23.4%的跨度内出现均匀分布的环形裂缝(图 6)，间隔为10 m左右.地表变形特征表明地震波的传播方向斜交于隧道轴线方向，容易引起隧道轴向拉伸、压缩变形，会导致横/斜/环向裂缝；加之致密安山岩和破碎安山岩成层交替分布，导致隧道呈现规律性环向裂缝，如图 14所示.
5 结论
1) 基于高精度航空机载激光测量数据(LiDAR)，通过分类迭代法技术获取受2016年日本熊本地震影响作用下俵山隧道周边地表变形分布特征.计算结果与日本国土地理院公布的实测值一致，表明该研究方法的有效性和准确性.
2) 研究分析了隧道周边地表空间变形和隧道破坏的关系，发现隧道周边的地表水平变形和隧道的破坏模式有较强的相关性，而垂直变形的相关性不大.
3) 强烈的地表水平变形运动特征能够在一定程度上反映地下结构震害的真实情况，如衬砌、路面和施工缝等各种结构的破坏特征.地表水平变形运动的方向能够在一定条件下(如排除地震断层等特定地质构造的影响等条件)反映地震波的传播方向.这一发现也验证了俵山隧道出现的均匀分布环向裂缝关于地震波传播方向的猜测.
4) 在一定程度上，掌握地震变形场的详细信息，特别是地表变形的方向，便于划分地下结构的震害影响区域，对于了解地下结构的地震力学反应和抗震减震具有一定的应用价值.
| [1] |
杜修力, 李洋, 许成顺, 等. 1995年日本阪神地震大开地铁车站震害原因及成灾机理分析研究进展[J]. 岩土工程学报, 2018, 40(2): 223. DU Xiuli, LI Yang, XU Chengshun, et al. Review on damage causes and disaster mechanism of Daikai subway station during 1995 Osaka-Kobe earthquake[J]. Chinese Journal of Geotechnical Engineering, 2018, 40(2): 223. DOI:10.11779/CJGE201802002 |
| [2] |
崔光耀, 王明年, 于丽, 等. 汶川地震公路隧道洞口结构震害分析及震害机理研究[J]. 岩土工程学报, 2013, 35(6): 1084. CUI Guangyao, WANG Mingnian, YU Li, et al. Seismic damage and mechanism of portal structure of highway tunnels in Wenchuan earthquake[J]. Chinese Journal of Geotechnical Engineering, 2013, 35(6): 1084. |
| [3] |
ZHANG X P, JIANG Y J, SUGIMOTO S. Seismic damage assessment of mountain tunnel: A case study on the Tawarayama tunnel due to the 2016 Kumamoto earthquake[J]. Tunnelling and Underground Space Technology, 2018, 71: 138. DOI:10.1016/j.tust.2017.07.019 |
| [4] |
AYDAN O. Rock dynamics[M]. CRC Press: Balkema, Netherlands, 2017.
|
| [5] |
YAMAZAKI F, LIU W. Remote sensing technologies for post-earthquake damage assessment: A case study on the 2016 Kumamoto earthquake[C]// 6th Asia Conference on Earthquake Engineering. Cebu City, Philippines, 2016
|
| [6] |
KAWAMOTO S, HIYAMA Y, OHTA Y, et al. First result from the GEONET real-time analysis system (REGARD): The case of the 2016 Kumamoto earthquakes[J]. Earth, Planets and Space, 2016, 68(1): 190. DOI:10.1186/s40623-016-0564-4 |
| [7] |
KAWAMOTO S, HIYAMA Y, KAI R, et al. Crustal deformation caused by the 2016 Kumamoto earthquake revealed by GEONET[J]. Bulletin of Geospatial Information Authority of Japan, 2016, 64: 27. |
| [8] |
MOYA L, YAMAZAKI F, LIU W. Calculation of coseismic displacement from lidar data in the 2016 Kumamoto, Japan, earthquake[J]. Natural Hazards and Earth System Sciences, 2017, 17(1): 143. DOI:10.5194/nhess-2016-315 |
| [9] |
STURZENEGGER M, STEAD D. Close-range terrestrial digital photogrammetry and terrestrial laser scanning for discontinuity characterization on rock cuts[J]. Engineering Geology, 2009, 106(3/4): 163. DOI:10.1016/j.enggeo.2009.03.004 |
| [10] |
SLOB S. Automated rock mass characterization using 3-D terrestrial laser scanning[D]. TU Delft: Delft University of Technology, 2010
|
| [11] |
GIGLI G, CASAGLI N. Semi-automatic extraction of rock mass structural data from high resolution LIDAR point clouds[J]. International Journal of Rock Mechanics and Mining Sciences, 2011, 48(2): 187. DOI:10.1016/j.ijrmms.2010.11.009 |
| [12] |
RIQUELMEA J, ABELLáN A, TOMáS R, et al. A new approach for semi-automatic rock mass joints recognition from 3D point clouds[J]. Computers and Geosciences, 2014, 68: 38. DOI:10.1016/j.cageo.2014.03.014 |
| [13] |
RIQUELMEA J, ABELLáN A, TOMáS R. Discontinuity spacing analysis in rock masses using 3D point clouds[J]. Engineering Geology, 2015, 195: 185. DOI:10.1016/j.enggeo.2015.06.009 |
| [14] |
DOWDING C H, ROZAN A. Damage to rock tunnels from earthquake shaking[J]. ASCE Journal of the Geotechnical Engineering Division, 1978, 104(2): 175. |
| [15] |
TAKAI S, DATE H, KANAI S, et al. Accurate registration of MMS point clouds of urban areas using trajectory[C]//ISPRS Annals of the Photogrammetry, Remote Sensing and Spatial Information Sciences, 2013, Ⅱ-5/W2: 277. DOI: 10.5194/isprsannals-Ⅱ-5-W2-277-2013
|
| [16] |
DEMANTKE J, MALLET C, DAVID N, et al. Dimensionality based scale selection in 3D LiDAR Point Cloud[J]. ISPRS Workshop Laser Scanning, 2011, 38(Part 5): W12. DOI:10.5194/isprsarchives-XXXVIII-5-W12-97-2011 |
| [17] |
ODA K, HATTORI S, TAKAYAMA T. Detection of slope movement by comparing point clouds created by SFM software[J]. ISPRS-International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, 2016, XLI-B5: 553. DOI:10.5194/isprsarchives-XLI-B5-553-2016 |
| [18] |
MUKOYAMA S, SATO T, TAKAMI T. Estimation of coseismic surface displacement in the Aso Caldera area before and after the 2016 Kumamoto earthquake by topographical data analysis from differential LiDAR DEM analysis[R]. Report of the 2016 Kumamoto, Oita Earthquake Disaster Research Mission. Japan Society of Engineering Geology, 2017: 55
|