地球重力场研究是间接研究构造应力场的有效方法^[1-6]。国内学者研究建立重力异常及垂线偏差转换水平构造应力场的方法^[2]，计算出不同区域的水平构造应力场^[2-4]；提出动态重力场转换动态水平构造应力场的理论方法^[5]；运用有限元求解华北地区重力势能差产生的水平构造应力^[6]。研究表明，通过地应力、静岩应力、构造应力的内在联系可以实现地球重力场数据向应力场数据的转换^[3]。包含龙门山断裂带、鲜水河-安宁河断裂带等区域的川西地区，是研究青藏高原向东南扩展变形及其与华南地块相互作用的关键地带，也是探讨汶川M_S8.0地震和芦山M_S7.0地震构造变形机制、开展地震风险评估的关键地域^[7-8]。本文选取构造活动强烈的川西地区重力观测数据，研究探讨水平构造应力场的转换。通过对川西地区(26.5°~32°N，100°~105°E)的静态构造应力场进行分析，探究构造应力场背景及构造运动与地震孕育之间的联系。

1 重力场与构造应力场的转换

应用希尔伯特三维位场转换公式^[2]，由重力异常g_z解出重力水平分量g_x和g_y：

$ \left\{ \begin{array}{l} {g_x} = \\ - \frac{1}{{2\pi }}\int_{ - \infty }^{ + \infty } {\int_{ - \infty }^{ + \infty } {\frac{{(x - \xi ){g_z}}}{{{{[{{(x - \xi )}^2} + {{(y - \eta )}^2} + {H^2}]}^{3/2}}}}{\rm{d}}\xi {\rm{d}}\eta } } \\ {g_y} = \\ - \frac{1}{{2\pi }}\int_{ - \infty }^{ + \infty } {\int_{ - \infty }^{ + \infty } {\frac{{(y - \eta ){g_z}}}{{{{[{{(x - \xi )}^2} + {{(y - \eta )}^2} + {H^2}]}^{3/2}}}}{\rm{d}}\xi {\rm{d}}\eta } } \end{array} \right. $

(1)

式中，(x，y)为计算点坐标，(ζ，η)为整个区域的流动坐标，H为空间延拓高度。

由式(1)得出重力总水平梯度的大小与方向：

$ \alpha = \arctan \frac{{{g_y}}}{{{g_x}}} $

(2)

$ {g_L} = {({g_x}^2 + {g_y}^2)^{1/2}} $

(3)

其中，重力总水平梯度的方向近似于水平最大主应力方向。

水平构造应力计算公式为^[2]：

$ \Delta {\sigma _{{\rm{xx}}}} = \frac{{g{\rho _{x, y}}{g_L}}}{{4\pi f{\rho _m}}} $

(4)

式中，ρ_m为地幔密度，取值3 300 kg/m³；f为万有引力常数，取值6.67×10^-11N·m²/kg²；g为重力加速度，取值9.8 m/s²；ρ_{x, y}为经过均衡改正的单位柱体密度：

$ {\rho _{x, y}} = {\rho _c} - {\sigma _{x, y}} + {(\sigma )_{x, y}} $

(5)

δ_{x, y}为地形补偿密度：

$ {\sigma _{x, y}} = \frac{{{\rho _c}{H_{x, y}}}}{D} $

(6)

(δ)_{x, y}为山根补偿密度：

$ {(\sigma )_{x, y}} = \frac{{({\rho _m} - {\rho _c})({d_{x, y}} - D)}}{D} $

(7)

式中，ρ_c为地壳平均密度，取值2 670 kg/m³；H_{x, y}为单位柱体地形的平均海拔高；d_{x, y}为单位柱体的山根平均深度；D为补偿面深度，Airy模型取30 km。

2 川西地区构造应力场计算 2.1 数据来源

本文利用EGM2008地球重力场模型^[9]数据计算静态构造应力场。该模型由地面区域平均重力数据、卫星测高数据、GRACE卫星重力数据等重力资料与全球高分辨率地形数据构成，可提供2 159次、2 190阶的球谐系数，空间分辨率为5′，即9 km分辨率。本文所使用的EGM2008地球重力场模型来源于斯克里普斯海洋研究所(SIO)2016年发布的V24.1重力模型和2014年发布的V18.1地形模型，数据空间分辨率为1′。其中，SIO V24.1重力模型由EGM2008的陆地重力异常数据及卫星测高的海域重力异常数据构成，新版EGM2008地球重力场模型对大陆重力场改动较小，主要对海洋重力场进行了优化^[10]，能较好地反映我国大陆因地形起伏引起的重力场变化^[11]；SIO V18.1地形模型由ETOPO1陆地地形数据与海洋测深及卫星测高所得到的海底地形数据构成^[10]，利用SIO V18.1地形模型绘制的川西地区地形图如图 1所示。

2.2 川西地区地壳浅层水平构造应力场

静态构造应力场计算时，积分区域经离散化、格网化后以2°×2°范围为积分区域进行计算。数据采用Kriging法插值格网。在重力场与构造应力场转换时，将选取范围外延2°，以消除边界效应。通过式(3)计算川西地区重力总水平梯度，利用式(4)计算川西地区地壳浅层水平构造应力场，绘制出川西地区地壳浅层水平构造应力矢量图(图 2)和地壳浅层水平构造应力等值线图(图 3)。图 2中箭矢长度表示水平构造应力大小，箭矢所指方向为箭矢起点处最大水平主应力的迹线方向。箭矢颜色越趋于红色表示水平构造应力越大，越趋于蓝色表示水平构造应力越小。红色折线表示活动断层，黑色折线表示块体分界线^[12]。图 3中越趋于红色表示水平构造应力场等值线越密集，构造应力值越大，越趋于蓝色表示水平构造应力场等值线越稀疏，构造应力值越小。红色五角星标示汶川M_S8.0地震和芦山M_S7.0地震震中位置。川西地区水平构造应力在0~17 MPa范围内。

2.3 川西地区构造应力场随深部变化结果

根据重力延拓法计算川西地区10 km、20 km和50 km深度处的重力总水平梯度，再由式(4)计算出不同深度处的水平构造应力场如图 4、图 5和图 6所示。蓝色五角星所在位置为相应深度水平构造应力最大值处，黑色菱形为汶川M_S8.0地震和芦山M_S7.0地震震中位置。

3 结果分析与讨论

由图 2可见，川西地区构造应力比较复杂，北部龙门山地区的构造应力总体呈东南方向，东南部的构造应力总体呈近东西向。在龙门山断裂带区域(30°~32°N，102.5°~105°E)，构造应力较大且聚集明显，表现为巴颜喀拉块体向东南方向运动时受到华南地块的强烈阻碍。在鲜水河-安宁河断裂带区域(29°~30°N，101.5°~102.5°E)，构造应力主要表现为近南北向，表明该地区受到北部块体的挤压。

此外，龙门山应力区目前处于近东西向的挤压构造应力环境中。在横向上，水平构造应力由巴颜喀拉块体的北东-南西向到龙门山块体的近东西向，再到华南地块的北西西-南东东向，表现为顺时针方向旋转；在纵向上，水平构造应力由华南地块的北西西-南东东向，到川滇块体的近东-西向，亦表现为顺时针方向旋转。应力方向在不同块体间的转变，反映出青藏高原东缘巴颜喀拉块体在南东向水平运动时对龙门山块体产生了挤压，同时由于华南地块特别是四川盆地的强烈阻挡，导致龙门山推覆构造带出现明显的地壳缩短和山体隆升。该结果与通过震源机制解获得的最大水平主应力优势方向较为一致^[13]。

由图 3可见，川西地区水平构造应力变化范围为0.050~16.958 MPa，这一结果与郭飞霄等^[4]利用重力数据反演川滇地区构造应力场的结果相一致。应力集中或构造运动活跃的区域，等值线密集程度越高。龙门山断裂带(30°~32°N，102.5°~105°E)、鲜水河-安宁河断裂带(29°~30°N，101.5°~102.5°E)等构造强活动区，对应于等值线图中红色区域；等值线稀疏或无等值线出现的地带为非构造应力聚集区，表明相应地带的构造运动不活跃，反映了华南地块尤其是四川盆地(28.5°~31.5°N，103°~105°E)为构造运动较稳定区域；等值线图中马尔康、道孚、峨边、雅江等地区附近的等值线密集程度介于二者之间，为构造运动相对活跃区。与孟文等^[13]获得的龙门山断裂带构造应力实测数据1.74~25.69 MPa相比，本文所得结果略小，可能是由于本文的构造应力仅来源于重力场数据反演，虽然体现了实际构造应力的大部分信息，但仍有部分构造应力信息缺失。但是，构造应力分布与实测数据较为一致，构造应力最大值均出现在康定附近，龙门山南段在较大尺度范围内地壳浅层构造应力强度较高，龙门山中段构造应力强度有所降低。

从图 2和图 3可以看出，川西地区恰好处于地壳活跃地带，历史上也曾发生过多次强烈地震，如汶川M_S8.0地震和芦山M_S7.0地震。由图可见，龙门山断裂带构造应力明显聚集，且构造应力数值高于其他地区，地壳浅层构造应力最大值接近17 MPa，说明该地区构造活动强烈。龙门山应力聚集区目前正处于近东西向的挤压构造应力环境之中，与该区域逆冲断层为主的强烈构造活动情况相一致。GPS速度场资料显示，汶川地震以前，龙门山块体运动速度由藏东块体到龙门山块体衰减了49%^[14]。汶川地震以后，巴颜喀拉地块与川滇块体相对于华南地块的运动速率进一步降低^[15]，而该地区的水平构造应力相对于周边地区呈急剧增强态势，说明此区域断层在漫长时期(至少为万年尺度以上)是强活跃的，目前出现强闭锁现象。鲜水河-安宁河断裂带与龙门山断裂带交汇处，尤其在康定南部有较大的应力聚集，构造应力方向主要沿断裂带方向分布，构造活动强烈，且主要为走滑断层。川西地区南部有一定程度的构造应力聚集，构造活动较明显。华南地块尤其是四川盆地，应力分布较平缓，应力聚集程度较低，反映出该地区构造运动相对稳定。

如图 4、图 5和图 6所示，各深度处构造应力最大值分别为：地壳浅层16.958 MPa、10 km深度11.388 MPa、20 km深度8.108 MPa、50 km深度4.484 MPa，构造应力大小随着深度增加呈递减趋势, 该结果与通过重力场算得的构造应力分布结果一致^[3-4]。此外，10 km、20 km和50 km深度的水平构造应力最大值所在位置均处于汶川M_S8.0地震和芦山M_S7.0地震震中附近，说明该地区不同深度构造活动均比较活跃，汶川M_S8.0地震和芦山M_S7.0地震在一定程度上均印证了这一特征。

图 7为2009-01~2017-09不同深度处地震统计图，信息来源于中国地震台网中心地震数据管理与服务系统川西地区4级以上地震的统计资料。综合图 4~图 7可得，随着深度的增加构造应力逐渐变小，构造运动相对稳定，地震活动频率降低。

图 8给出了川西地区地壳浅层构造应力场与地震分布情况，信息来源于中国地震信息网1900~2017年川西地区6级以上地震的统计资料。图中红色圆点标示地震分布，圆点越大表示震级越大。从图 8可以看出，在构造应力场等值线相对密集的地区构造运动活跃、地震活动频发，表明地震活动与构造应力聚集分布存在一定的对应关系。龙门山断裂带与鲜水河-安宁河断裂带属水平构造应力较大的地区，中强、强地震分布较为稀疏，可能会导致断层闭锁而造成能量积累，仍有强震发生的可能。需要注意的是，图中标示的地震几乎均发生于应力突变区域，可见应力分布不均匀以及应力突变可能是导致地震频发的一个重要因素。

4 结语

1) 利用重力异常数据能够快速获取水平构造应力场。获得的水平构造应力场能够基本反映川西地区地壳构造活动状况，川西地区构造活动特征可通过水平构造应力大小及方向较好地诠释。然而，通过重力资料获得的水平构造应力场虽能反映大部分水平构造应力特征，但仍有部分水平构造应力信息缺失。

2) 以静态水平构造应力场作为背景构造应力场，结合流动重力观测资料、GPS观测位移资料、地震波观测资料等，进一步开展动态水平构造应力场研究，能够更深入地探究构造应力场及重力场与地震孕育的关系，为研究断层运动机理、探索孕震过程作好铺垫。