应用非连续模型研究断层对地壳应力的影响

雷显权，陈运平

(中南大学 地球科学与信息物理学院 计算地学研究中心，湖南 长沙，410083)

摘要：通过有限元数值模拟的方法，应用非连续的接触模型研究断层对地壳应力的影响。研究结果表明：非连续模型比连续模型能够更好地模拟断层，突出了断层对地壳应力的影响作用，避免了单纯使用软弱带模拟断层所导致的应力削弱现象；非连续模型通过接触单元可以描述连续模型所不能反映的断层带与围岩之间的不连续错动；断层端部和复合断层交叉处附近是应力集中区；断层带具有较高的剪切应力，但其自身及附近局部区域的最大主应力受到了一定的弱化；断层的弹性模量、断层与区域主应力的夹角以及交叉断层是地壳应力的重要影响因素，平行断层是弱影响因素，断层泊松比对地壳应力的影响不大。

关键词：非连续模型；断层；地壳应力；数值模拟

中图分类号：P553          文献标志码：A         文章编号：1672-7207(2011)08-2379-08

Effects of faults on crustal stress using discontinuous model

LEI Xian-quan, CHEN Yun-ping

(Research Centre of Computational Geosciences, School of Geosciences and Info-Physics,

Central South University, Changsha 410083, China)

Abstract: The effects of faults on crustal stress applying discontinuous contact model were studied through FEM numerical simulation. The results show that discontinuous model is better than continuous model in modeling faults, and it highlights the influencing effect of faults on crustal stress and improves stress-weaking resulting from simply modeling faults with weak-belts. The discontinuous dislocation between fault belt and surrounding rock can be depicted by discontinuous model rather than continuous model. The ends of faults and the vicinity of the cross of complex faults are stress concentrational region. The shear stress in fault belt is high but the maximum principal stress near faults in itself is weakened to some degrees. The elastic modulus of faults, the angle between faults strike and regional maximum principal stress and intersectant faults are important influencing factors for crustal stress. Parallel faults are weak factor. The poisson ratio has little influence on crustal stress.

Key words: discontinuous model; faults; crustal stress; numerical modeling

地壳应力是地震预测和地震危险性评价的一个重要依据，也是地质环境、地质工程、岩土工程设计和施工的重要基础参数^[1-2]。地质构造运动会引起地壳形变和应力，地应力局部集中会造成地震活动和岩体破坏。若岩体内存在断裂等不连续面, 则其附近的应力状态会发生复杂的变化^[3-5]。地应力是一个地区的区域稳定性和岩体稳定性的重要因素，它对建筑物、构筑物和地下工程的安全有直接的影响。围岩发生失稳破坏大多沿岩体不连续面即节理、断层或者附近进行。在岩体高应力区，边坡开挖、隧洞穿凿、巷道采矿、地下空间工程等岩体开挖会引起一系列与应力释放相关的岩体变形和破坏现象，恶化岩体的工程地质条件，直接威胁到施工工程的安全和稳定^[6-8]。因此，对地应力的研究受到工程地质界的广泛重视。在工程力学模拟计算中，往往采用连续变形算法，这对断裂构造的模拟有很大的不足，在理论计算和工程分析中都难免出现相当大的误差甚至错误^[9]，特别是当存在一些大断裂和断裂两侧地壳变形跃变剧烈时，仍然采用连续变形算法是不够的。计算机技术和计算方法的发展，不连续体概念和不连续变形理论逐渐成熟，为工程力学数值模拟和计算地球科学提供了新的方法^[10-13]。在此，本文作者采用非连续接触单元模拟断层，根据库仑摩擦准则判断剪切滑动和计算剪切应力，通过有限元数值模拟方法来探讨断层各要素对地壳应力的   影响。

1  理论和方程

地壳岩石往往被断层分割成众多离散的块体，块体内部的变形可以采用固体力学方程来描述，而不同块体彼此之间则表现为断层接触关系。模型遵循的物理力学方程如下^[14]。

(1) 平衡方程。物体平衡的力学条件为：

             (1)

(2) 几何方程。物体在应力作用下产生形变，应变和位移满足如下几何关系：

           (2)

(3) 本构方程。对于各向同性线弹性材料，应力与应变对载荷的响应符合胡克定律：

           (3)

其中：σ_ij为应力；ε_ij为应变；f_i为体力；x_i为坐标；u_i为位移；θ为体应变；λ和μ为拉梅系数，λ=Eν/[(1+ν)(1-2ν)]，μ=E/[2(1+ν)]；E和ν分别为弹性模量和泊松比；δ_ij为Kronecker变量，若i=j，则δ_ij=1；若i≠j，则δ_ij=0。

(4) 接触应力方程。断层往往不只是一个简单的错动面，而是相对围岩稍微软弱的断层带。断层带和围岩自身的变形是连续的，但二者之间是不连续的接触关系。由于断层带和围岩的接触区域通常并不完全确定，所以，本文模型采用点/面接触方式来处理二者的接触关系。将断层带与围岩的接触边界相互视为对方的接触面和目标面，从而组成接触对，然后，通过接触力学方法分析断层的受力和运动状态^[15]。接触分析需要计算垂直于目标面的法向接触应力和平行于目标面的切向接触应力，采用Pinball算法进行计算。

接触面和目标面之间的间隙(穿透量)用g来表示，当接触面穿过目标面时，就发生接触穿透，规定此时g为负值。采用罚函数法计算法向接触应力：

              (4)

其中：f_n为法向接触应力；K_n为法向接触刚度。

切向接触应力是由接触面在目标面上移动所产生的摩擦力引起的，若接触面沿目标面的切向位移的弹性分量为，则切向接触应力为：

      (5)

其中：f_t为切向接触应力；K_t为切向接触刚度；F为静态/动态摩擦因子；τ为库仑滑动摩擦力：τ=κf_n；κ为摩擦因数。

2  模型和边界条件

2.1  基本假设

对于发育断层的二维平面地壳模型，进行以下假设：

(1) 断层带与围岩都为各向同性的线弹性材料；

(2) 断层带相对围岩较软弱；

(3) 断层带与围岩为非连续的接触关系；

(4) 断层带与围岩之间的相对滑动符合库仑摩擦定律；

(5) 顶板和底板的垂向刚度极大而水平剪切刚度极小。

2.2  模型的建立

采用二维平面模型，考虑水平面上的正方形岩体和其中的不同几何形态的非贯通型断层。岩体长×宽为1 km×1 km，断层带宽度为10 m。对于具有一定宽度的断层带，非连续模型将其外边界处理为与岩体紧挨的接触面。

2.3  模型力学参数

真实材料参数本应从实验中获取，但不同时代、不同地区、不同深度、不同性质岩石的参数都不尽相同。本文进行一般性讨论，岩体和断层带的密度、弹性模量、泊松比和摩擦因数在合理范围内取值^[16]。对于接触力学参数，如法向接触刚度因子、切向接触刚度因子和静态/动态摩擦因子，本文均取1.0，这适合于ANSYS求解大多数接触问题^[15]。岩体和断层接触面的摩擦因数为0.4。除了在分析断层带力学参数对地壳应力的影响时需要改变断层带的弹性模量和泊松比外，所有模型均采用表1所示的参数。

表1  有限元模型力学参数

Table 1  Mechanical parameters of finite element model

2.4  边界条件

模型受到区域最大主应力σ₁和最小主应力σ₃的挤压作用，σ₁=10 MPa，σ₃=5 MPa。模型左边界和前部边界在法向上被约束，右边界施加最大主应力σ₁，后部边界施加最小主应力σ₃。所有边界剪切自由。

岩体中发育的1条平直断层模型的有限元分网及边界条件如图1所示。其他模型有与此相同的边界条件，只是断层几何形态有所不同。

图1  有限元模型及边界条件

Fig.1  Finite element model with boundary conditions

3  计算结果及分析

3.1  连续模型和非连续模型地壳应力状态的比较

连续模型与非连续模型的区别在于后者考虑了断层带与围岩之间不连续运动的接触关系，在结合处通过接触单元把二者联系起来，二者在接触位置可以发生错动；而对于连续模型，断层带与围岩在结合处具有公共节点，因此，二者之间不会发生错动(图2)。

图2  连续模型与非连续模型示意图

Fig.2  Schematic diagrams of continuous and discontinuous model

采用相同的模型参数和边界条件分别计算连续模型和非连续模型的地壳应力，结果表明无论是连续模型还是非连续模型，断层及围岩的应力状态表现出3个相同的特点：(1) 弱化了断层带两端局部较小区域的最大主应力；(2) 断层带具有较高的剪切应力；    (3) 最大主应力和剪切应力的高值都位于断层端部。但是，非连续模型在更大程度上明显地影响了断层端部的最大主应力和剪切应力分布状态。沿着断层的最大主应力和剪切应力进行对比，发现非连续模型比连续模型更加突出了断层对其附近地壳应力的影响作用(图3)。上述特征表明：非连续模型能够更好地模拟断层，有效克服了单纯使用软弱带模拟断层而不考虑其不连续运动所导致的应力削弱的问题。

3.2  连续模型和非连续模型断层走滑变形特征的  比较

断层走滑量可以通过计算断层两侧与断层带接触的岩体单元节点在断层走滑方向上的位移之差得到，而且可以根据位移关系推测断层的走滑性质。图4所示为连续模型和非连续模型断层走滑方向上的位移沿1条垂直于断层的跨断层测线的变化(测线见图1，方向为左上—右下)。从图4可以看出：在给定的区域构造应力作用下，连续模型和非连续模型都反映出断层发生了左旋走滑作用，但二者给出的走滑量相差很大：非连续模型断层走滑量为1.73 cm，而连续模型的走滑量仅为0.32 cm；非连续模型断层走滑变形由3段构成，其中第1和第3段为断层带与围岩之间的不连续错动，第2段为断层带内部的连续变形，且走滑量主要为不连续错动所贡献；而连续模型断层仅有断层带的连续变形这一部分，由于断层带与围岩在结合处具有公共节点，故二者并没有发生错动。野外地质观察到的大量擦痕和阶步是断层不连续错动的具体体现，而非连续模型可以通过接触单元描述断层带与围岩之间的这种不连续变形，进一步表明非连续模型比连续模型能够更好地模拟断层。

图3  连续模型和非连续模型地壳应力沿断层的变化

Fig.3  Changes of crustal stress along fault of continuous and discontinuous model

图4  连续模型和非连续模型断层走滑方向上的位移沿跨断层测线的变化(测线见图1)

Fig.4  Changes of displacement in fault strike-slipping direction along test line A-B striding over fault

of continuous and discontinuous kinematic model (test line shown in Fig.1)

由于断层端部是地壳应力的高值分布区即应力集中区，该处应力状态能够较好地反映断层对地壳应力的影响。

3.3  断层力学参数对地壳应力的影响

3.3.1  弹性模量的影响

采用相同的非连续模型计算只改变断层弹性模量时地壳应力的变化，结果如图5所示。从图5可以看出：随着弹性模量的增大，断层端部的最大主应力和剪切应力明显减小，尽管最大主应力并没有像剪切应力那样一直递减，但数值的变化是显而易见的；断层中部的最大主应力变化很大，但剪切应力的变化较小；断层远处的最大主应力和剪切应力的变化都很小，表明断层端部地壳应力对断层的弹性模量非常敏感。确切地说，单独改变断层弹性模量实际上是改变了围岩与断层弹性模量的差别。当断层弹性模量取较小值时，断层与围岩的弹性模量的差别较大，断层端部具有较大的最大主应力和剪切应力。因此推断：在相同的区域主应力作用下，围岩与断层的弹性模量相差越大，断层端部的地壳应力就越大。

3.3.2  泊松比的影响

采用相同的非连续模型计算只改变断层泊松比时地壳应力的变化，如图6所示。从图6可以看出：断层中部和远处的最大主应力和剪切应力几乎都不随泊松比的改变而变化。尽管断层端部应力有随泊松比增大而减小的微弱趋势，但变化量非常小，表明断层泊松比对地壳应力的影响是微弱的。

图5  非连续模型断层附近地壳应力与断层弹性模量的关系

Fig.5  Relationship between crustal stress near faults and elastic modulus of faults for discontinuous model

图6  非连续模型断层附近地壳应力与断层泊松比的关系

Fig.6  Relationship between crustal stress near faults and poisson ratio of faults for discontinuous model

3.4  断层与区域主应力具有不同夹角时断层端部地壳应力的变化

采用相同的力学参数计算改变断层与区域最大主应力的夹角(θ)时断层端部地壳应力的变化，如图7所示。从图7可以看出：断层端部最大主应力和剪切应力都随θ的改变而变化。最大主应力在θ约为45°时具有极大值，剪切应力在θ约为25°和70°时具有极大值。尽管应力极大值对应的具体角度可能有误差，但在一定程度上反映了断层端部地壳应力与断层和区域主应力之间夹角的变化关系。连续模型与非连续模型有相似的变化特征，只是前者在数量上比后者的小。

3.5  复合断层对地壳应力的影响

在地质构造活动强烈的地区，复合断层是很常见的，根据断层的发育关系大致可归纳为2类：近似平行断层和交叉断层。平行断层和交叉断层对地壳应力的影响采用非连续模型进行计算和讨论。

3.5.1  平行断层的影响

考虑2条长、宽相同的平行断层，改变其中一条断层与另一条固定断层的距离，计算固定断层端部地壳应力的变化，如图8所示。从图8可以看出：最大主应力和剪切应力都随断层间距减小而逐渐增大，但是，当断层间距较大时，应力的变化非常小，此时平行断层对地壳应力的影响几乎可以忽略。仅在断层间距很近时，应力随距离减小而快速增大。但是，无论是最大主应力还是剪切应力，应力的绝对改变量都并不是很大，表明平行断层是地壳应力比较弱的影响因素。

图7  断层端部地壳应力与断层和区域最大主应力之间夹角θ的关系

Fig. 7  Relationship between crustal stress at end of faults and angle θ between fault-strike and regional maximum principal stress

图8  非连续模型断层端部地壳应力与平行断层间距d的关系

Fig.8  Relationship between crustal stress at end of faults and distance d between two parallel faults for discontinuous model

3.5.2  交叉断层的影响

交叉断层的宏观形态多为“X”形，计算“X”形交叉断层具有不同夹角(β)时的地壳应力，结果表明：除了断层端部，断层交叉处也是应力的集中区，这与孤立断层的地壳应力状态明显不同；应力扰动范围随交叉断层夹角的变化而发生改变；在断层交叉处附近，位于区域最大主应力(σ₁)方向上的2个三角区的最大主应力表现为应力强化，而位于区域最小主应力(σ₃)方向上的2个三角区的最大主应力表现为应力弱化；断层交叉处及附近的剪切应力表现为高值且呈中心对称形态。图9显示了地壳应力与交叉断层夹角(β)的关系。从图9可以看出：交叉断层夹角对地壳应力的影响  很大。

图9  非连续模型断层交叉处地壳应力与断层夹角β的关系

Fig.9  Relationship between crustal stress at faults crossing point and intersecting angle β of two faults for discontinuous model

4  结论

(1) 非连续模型比连续模型能够更好地模拟断层，突出了断层对地壳应力的影响作用，克服了单纯使用软弱带模拟断层所导致的应力削弱现象。

(2) 非连续模型通过接触单元可以描述连续模型所不能反映的断层带与围岩之间的不连续错动。

(3) 断层端部和复合断层交叉处附近是应力集中区；断层带具有较高的剪切应力，但其自身及附近局部区域的最大主应力受到了一定的弱化。

(4) 断层的弹性模量、断层与区域主应力的夹角和交叉断层是地壳应力的重要影响因素，平行断层是弱影响因素，断层泊松比对地壳应力的影响不大。

参考文献：

[1] 谢富仁, 陈群策, 崔效锋, 等. 中国大陆地壳应力环境基础数据库[J]. 地球物理学进展, 2007, 22(1): 131-136.
XIE Fu-ren, CHEN Qun-ce, CUI Xiao-feng, et al. Fundamental database of crustal stress environment in continental China[J]. Progress in Geophysics, 2007, 22(1): 131-136.

[2] 张咸恭, 王思敬, 张倬元. 中国工程地质学[M]. 北京: 科学出版社, 2000: 97-103.
ZHANG Xian-gong, WANG Si-jing, ZHANG Zhuo-yuan. Chinese engineering geology[M]. Beijing: Science Press, 2000: 97-103.

[3] 苏生瑞, 朱合华, 王士天, 等. 断裂物理力学性质对其附近地应力场的影响[J]. 西北大学学报: 自然科学版, 2002, 32(6): 655-658.
SU Sheng-rui, ZHU He-hua, WANG Shi-tian, et al. The effect of fracture properties on stress field in the vicinity of a fracture[J]. Journal of Northwest University: Natural Science Edition, 2002, 32(6): 655-658.

[4] 朱维申, 阮彦晟, 李晓静, 等. 断层附近应力分布的异常和对隧洞稳定的影响[J]. 地下空间与工程学报, 2008, 4(4): 685-689.
ZHU Wei-shen, RUAN Yan-sheng, LI Xiao-jing, et al. Abnormal stress distribution adjacent to a fault and it’s influence on stability of tunnel[J]. Chinese Journal of Underground Space and Engineering, 2008, 4(4): 685-689.

[5] 孙礼健, 朱元清, 杨光亮, 等. 断层端部及附近地应力场的数值模拟[J]. 大地测量与地球动力学, 2009, 29(2): 7-12.
SUN Li-jian, ZHU Yuan-qing, YANG Guang-liang, et al. Numerical simulation of ground stress field at ends and vicinity of a fault[J]. Journal of Geodesy and Geodynamics, 2009, 29(2): 7-12.

[6] 伍佑伦, 王元汉. 不连续面激活与巷道围岩破坏区关系的探讨[J]. 岩土力学, 2007, 28(6): 1197-1200.
WU You-lun, WANG Yuan-han. Discussion on relationship between activation of discontinuities and disturbed zone in surrounding rock mass of tunnels[J]. Rock and Soil Mechanics, 2007, 28(6): 1197-1200.

[7] 荣冠, 朱焕春, 王思敬. 锦屏一级水电站左岸边坡深部裂缝成因初探[J]. 岩石力学与工程学报, 2008, 27(S1): 2855-2863.
RONG Guan, ZHU Huan-chun, WANG Si-jing. Primary research on mechanism of deep fracture formation in left bank of Jinping first hydropower station[J]. Chinese Journal of Rock Mechanics and Engineering, 2008, 27(S1): 2855-2863.

[8] 晏长根, 刘彤, 伍法权. 复杂条件下大型地下洞室群的变形稳定性分析[J]. 工程地质学报, 2008, 16(1): 84-88.
YAN Chang-gen, LIU Tong, WU Fa-quan. Deformation and stability analysis of large-scale underground cavity group under complicated geological conditions[J]. Journal of Engineering Geology, 2008, 16(1): 84-88.

[9] 周维垣, 杨强. 岩石力学数值计算方法[M]. 北京: 中国电力出版社, 2005: 282-284.
ZHOU Wei-yuan, YANG Qiang. Numerical methods for rock mechanics[M]. Beijing: Chinese Electric Press, 2005: 282-284.

[10] 石根华. 块体系统不连续变形数值分析新方法[M]. 北京: 科学出版社, 1993: 50-80.
SHI Gen-hua. New numerical analysis methods for discontinuous deformation of block system[M]. Beijing: Science Press, 1993: 50-80.

[11] 党亚民, 陈俊勇, 晁定波. 中国及其邻区地壳应力场数值特征的研究[J]. 测绘科学, 2001, 26(2): 11-14.
DANG Ya-min, CHEN Jun-yong, CHAO Ding-bo. Numerical study on stress characteristic in China[J]. Science of Surveying and Mapping, 2001, 26(2): 11-14.

[12] 郑勇, 陈颙, 傅容珊, 等. 应用非连续性模型模拟断层活动对青藏高原应力应变场的影响[J]. 地球物理学报, 2007, 50(5): 1398-1408.
ZHENG Yong, CHEN Yong, FU Rong-shan, et al. Simulation of the effects of faults movement on stress and deformation fields of Tibetan Plateau by discontinuous movement models[J]. Chinese Journal of Geophysics, 2007, 50(5): 1398-1408.

[13] Bird P, Kong X H. Computer simulations of California tectonics confirm very low strength of major faults[J]. Geol Soc Am Bull, 1994, 106: 159-174.

[14] 傅容珊, 黄建华. 地球动力学[M]. 北京: 高等教育出版社, 2001: 80-90.
FU Rong-shan, HUANG Jian-hua. Geodynamics[M]. Beijing: Higher Education Press, 2001: 80-90.

[15] ANSYS Inc. Release 11.0 documentation for ANSYS: Contact technology guide[R]. Canonsburg PA:ANSYS Inc, 2007: 1-10.

[16] Turcotte D L, Schubert G. Geodynamics[M]. Cambridge: Cambridge University Press, 2002: 433-436.

(编辑 陈灿华)

收稿日期：2010-07-11；修回日期：2010-09-25

基金项目：高等学校博士学科点专项科研基金资助项目(200805331115)；国家自然科学基金资助项目(40874048)；湖南省研究生科研创新基金资助项目(CX2010B051)

通信作者：陈运平(1970-)，男，海南乐东人，副教授，从事地球动力学及计算地球科学的研究；电话：0731-88877077；E-mail：chyp@csu.edu.cn