有限差分数值模拟在跨断层资料中的应用

荆红亮，郭国祥，靳宝萍，张聪聪，杨 静，李 军，韩 磊

(1. 山西省地震局临汾地震监测中心站，山西 临汾 041000； 2.太原大陆裂谷动力学国家野外科学观测研究站，山西 太原 030025)

0 引言

跨断层测量是研究断层现今活动特征的重要手段之一，观测数据蕴含了区域和断层带地球动力学和地震孕育发展变化的宝贵信息[1-5]。在地表进行的跨断层形变观测受气象、地下水位、荷载变化等干扰影响[6]，此外还受核幔耦合(地球自转速率变化)作用的调制[7]，信息成分比较复杂。其中载荷变化对跨断层数据的干扰，是目前形变观测中面临的关键问题。近年来，随着城镇建设的快速发展，工程建设对观测场地的干扰越来越严重，这些干扰从根本上归结为载荷对场地的影响[8]。为获得完整、准确地反映断层运动的地震前兆信息,研究场地周围的荷载对近场的影响十分必要[9]。目前，有限差分FLAC 3D数值模拟方法，大量被运用到载荷对周边变形观测影响的理论分析中。如，李宏伟[10]对郭庄水库库水位对地倾斜场的模拟计算,朱航[11]在耿达短水准资料异常辨别中的应用等。

在异常核实的过程中，通过定量分析，能够判定异常变化是否为干扰所致。该文以临汾峪里跨断层Y5-Y4测段高差2018年6月至2020年9月出现的累计8.62 mm的张性变化，以及广胜寺测线2017年8月至2018年4月(G1-G2)高差曲线加速张性变化为研究主体，建立有限元数值模拟模型，根据场地地质构造情况，设置基本计算参数及边界条件，使用岩土力学工程专用数值计算程序FLAC 3D软件，对峪里房屋荷载对Y5点点位的影响，和广胜寺观测场地附近施工与G1点点位的变化进行定量估算,分析建筑物、施工与峪里、广胜寺测值变化之间的关系，为异常判定提供依据。

1 跨断层资料变化特征分析

第33页图1为峪里跨断层场地Y5-Y4测段(场地测线布设见第33页图2)1980年1月至2020年9月高差曲线图，曲线上升为张性运动。该测段未跨越断层，Y5点为土层点，自1980年观测以来，曲线呈继承性加速张性变化，速率变化平稳，在±1.0附近波动。2018年6月至2020年9月出现大幅张性变化，累计8.62 mm。工作人员及时进行异常核实，现场对观测仪器、标尺检验，复测结果均正常。对观测环境调查发现，Y5端点周围除东侧、西南侧为空地外，其余几侧为民房建筑，房屋建于2009-2010年，自2013年开始出现裂缝，近几年破损尤其严重，距Y5端点南侧15 m左右民房下沉比较严重，外墙和内部隔断墙都存在断裂缝，最大的有2～3 cm左右。在该背景下，分析认为高差异常变化是Y5土层点快速下降导致。

图1 峪里跨断层1980年至2020年高差曲线图Fig.1 Curve of height difference between 1980 and 2020 of Yuli cross fault

图2 峪里跨断层测线场地布设示意图Fig.2 Site layout of Yuli cross fault survey line

第33页图3为广胜寺跨断层场地G1-G2测段(场地测线布设见第33页图4)1980年1月至2020年9月高差曲线图。该测段跨越断层，G1点位为土层点。从多年观测资料分析，2011年9月至2015年4月高差曲线呈现缓慢压性变化，2015年5月至2017年2月曲线呈上下波动，速率平缓，2017年10月至2018年4月出现短期加速张性变化，变幅为1.21 mm。经过现场调查，观测场地内距离G1点位20 m左右有建筑施工的相关信息，且有压路机作业。该建筑工程从2017年11月开工,2018年1月结束,历时约3个月。

图3 广胜寺跨断层1980年至2020年高差曲线图Fig.3 Height difference curve of Guangsheng Temple crossing fault from 1980 to 2020

图4 广胜寺跨断层测线场地布设示意图Fig.4 Site layout of Guangsheng Temple cross fault survey line

2 荷载模型建立与计算

2.1 FLAC 3D 软件

岩土力学工程专用数值计算程序FLAC 3D软件，是在满足基本方程(平衡方程、几何方程、本构方程)和边界条件下推算建立的。有限差分法的基本思想是把连续的定解区域用有限个离散点构成的网格来代替, 把连续定解区域上的连续变量的函数用在网格上定义的离散变量函数来近似；求解有限差分方程组就可以得到原问题在离散点上的近似解, 然后利用插值方法从离散解得到定解问题在整个区域上的近似值解[12]。该程序能较好地模拟计算三维岩、土体及其他介质中工程结构的受力与变形形态，特别是岩、土体的动力稳定性分析、土与结构的相互作用分析，获得计算过程中节点、单元的参数，如坐标、位移、应力、应变和不平衡力等[13]。该文利用有限差分计算程序 FLAC 3D，结合峪里、广胜寺两个观测场地现场实测资料，建立有限差分法计算模型，根据计算模型所处的地形和荷载条件,只考虑在荷载作用下沉降变形及应力影响深度,指定特殊的边界条件,自动进行参数分析，计算出两个观测场地不同类型载荷对跨断层数据的位移量，分析与实际测量结果是否吻合。如果数值模拟结果同跨断层实际观测曲线发展趋势一致，说明观测数据受到荷载影响；如果相差较大，排除异常观测数据受荷载干扰，该方法可以应用于异常核实工作中。

2.2 峪里模型参数及计算

根据对周边地质调查的结果，峪里测线Y5点位周边的岩土体为一层，主要是砂土，4 m以下为粘性砂土。第34页图5中①②③为荷载区域，为简化计算，主要根据图5中实际测量结果，测点周围以出现裂缝的房屋为荷载，建立三维的有限元摩尔-库仑计算静力荷载模型。由于计算区域为不均匀型，需建立较密的网格，设置模型边界条件长为102 m，宽173 m，限制侧面水平位移。该分析的目的主要是为得到周边受约束的弹性材料，在上表面受均匀荷载时所表现出来的力学响应特征。本模型用fix 命令对边界节点的X、Y、Z 三个方向的速度进行约束，并考虑岩土体由自身重力产生的沉降，运用solve 命令对求解结果保存文件，同时运用plot 命令输出Z 方向的沉降量云图。观测场地岩土力学参数[14]如表1所示。

表1 峪里岩土体材料力学参数Table 1 Material mechanical parameters of Yuli rock and soil mass

图5 峪里基础计算资料图Fig.5 Calculation data chart of Yuli Foundation

如图6所示，通过模型计算得到房屋荷载最大沉降量为9.758 mm，利用pri命令得出有裂缝房屋自重荷载引起Y5点位发生的竖向位移量为1.296 mm，二者的量值存在较大的差异，得出，房屋荷载自重沉降并不是引起点位下沉的主要因素。跨断层测量的实际观测值为8.62 mm，此与荷载本身自重发生的沉降量非常接近。分析认为沉降的原因可能与房屋、Y5点位所处区域的地形有关。调查发现，存在裂缝的几处房屋以及Y5点位刚好处于边坡坡后的拐角处，是否为边坡失稳导致坡体发生变形而产生的沉降，有待于进一步的研究。

图6 峪里Y5点计算结果沉降量云图Fig.6 Cloud chart of settlement amount calculated by Y5 point in Yuli

2.3 广胜寺模型参数及计算

为确定广胜寺G1点位是否受到测区施工压路机作业的影响，建立三维的有限元摩尔-库仑计算动力荷载模型。假设测点底部固定，顶部受到周期荷载的作用，以此来分析动力荷载作用下的测点桩体单元的动力响应。图7为广胜寺G1测点附近基础计算资料图，为反映动力响应特征，设置模型边界为50 m×50 m，测桩单元长度4 m，选择地面18 t压路机，岩土体为砂土和黏土两层，岩土体材料力学参数如表1。利用函数wave=sin(pi*dytime/16.0)[15]进行数值模拟(Pi为常数π，dytime为时间)。图8为广胜寺G1测点沉降量云图，根据荷载理论模型计算得出的结果为0.001 mm，远小于实际下降1.21 mm的形变量，因此认为测区施工压路机作业，不会引起广胜寺测段G1-G2高差曲线短期变化，排除测区施工对测点的影响，可能是存在其他因素或者构造运动引起了地面形变。

图7 广胜寺G1测点附近基础计算资料图Fig.7 Data chart of foundation calculation near G1 measuring point of Guangsheng Temple

图8 广胜寺G1点计算结果沉降量云图Fig.8 Cloud chart of settlement of G1 point in Guangsheng Temple

3 结论与讨论

(1) 通过有限元差分模拟计算峪里、广胜寺观测场地内荷载及施工对测点点位的数值解，结果显示，峪里房屋荷载对测点的影响值为1.296 mm的沉降，与实际观测值相差较远，因此判定Y5点位下沉，不是房屋荷载自重沉降引起；广胜寺观测场地内压路机施工引起的沉降量为0.001 mm，远小于实际1.21 mm的下降形变量，因此G1-G2测段高差变化与场地内施工无关。

(2) 利用FLAC 3D 有限差分法模拟荷载干扰以及计算点位的沉降量，该方法比较直观、简便，能够定量分析荷载对地形变的影响，在日常出现变化而进行异常核实研究中具有一定的意义。