2015, Vol. 58
紫坪铺水库位于岷江上游，于2001开始动工修建，2004年底开始蓄水.2008年汶川地震后，因其距离汶川地震的震中不足10 km，且水库库尾在主发震构造带上，引起了人们的注意，是否水库蓄水触发了大地震的发生？雷兴林等(2008)通过解析方法计算求解紫坪铺水库蓄水载荷和渗流作用，认为紫坪铺水库对龙门山中央断层和山前断层有明显的作用，在震源深度处库仑应力变化可达几个kPa，满足触发汶川地震发生的条件.2009年1月，Science刊物在News of The Week 栏目中发表了曾获行星科学新闻奖的Kerr和Stone(2009)的“A Human Trigger for the Great Quake of Sichuan?”的报道评述，报道了紫坪铺水库可能诱发了汶川地震.之后，各个媒体也跟进，甚至报道是水库制造了杀人的地震，汶川地震是紫坪铺水库蓄水引起的(http:// www.telegraph.co.uk/news/worldnews/asia/china/ 4434400/Chinese-earthquake-may-have-been-man- made-say-scientists.html)等. 之后，陈颙(2009)定性地讨论了汶川大地震不符合水库地震发生的一般特征，截然反对汶川地震诱发大地震的观点.随后，众多研究学者进一步对紫坪铺水库是否诱发了汶川地震的发生进行研究分析.陈厚群等(2008)认为紫坪铺水库蓄水对北川—映秀断裂原有的水文地质条件没有产生影响及水库蓄水前后地震活动性与库水位不存在相关关系，紫坪铺水库放水时间与汶川地震的发生仅仅是一种巧合.马文涛等(2011)据紫坪铺水库地震台资料，应用双差法方法研究了汶川地震前震、主震和余震，得出紫坪铺水库与汶川地震的发生有密切关联.程万正等(2010)对紫坪铺水库地震台网和地震台记录的地震资料进行分析，认为水库库区的小震活动可以视为诱发的小震，但仍属于构造本身长期地震活动的一种起伏，汶川地震为巨大构造地震.
随着争论的深入，越来越多的研究者进行了定量的计算模拟和讨论.Ge等(2009)建立了2D均匀介质模型，采用数值模拟方法，分别计算了紫坪铺水库水体静态荷载和孔隙压力扩散情况下的库仑应力变化(-0.01~0.05 MPa)，据其计算结果得出紫坪铺水库蓄水使得汶川地震提前10~100年到来.张贝和石耀霖(2010)应用3D线弹性Boussinesq解析解计算了紫坪铺水库水载荷作用，蓄水使得逆掩断层稳定，而放水增加了危险性，紫坪铺水库在震前放水(从875~817 m)会引起汶川地震震源处的库仑应力变化为0.001 MPa.Deng等(2010)通过3D解析解和2D数值模拟计算，得出紫坪铺水库蓄水对大地震的发生影响很小.周斌等(2010)基于实际库区地质岩体性质，讨论了在孔弹性介质下的2D有限元模拟结果.Zhu等(2011)应用LBM方法求得紫坪铺水库蓄水引起汶川地震震源处库仑应力变化为0.022 MPa，紫坪铺水库可能触发汶川地震的发生.Gahalaut等(2010)采用3D解析方法，分别计算了紫坪铺水库蓄水造成的弹性应力场和孔压扩散场，认为紫坪铺水库蓄水与汶川地震的发生没有关系.Klose(2012)求得在12 km处紫坪铺水库蓄水影响库仑应力变化约为0.004 MPa，认为水库触发了大地震的发生.孙玉军等(2012)利用3D孔隙弹性解耦模型探讨紫坪铺水库对汶川地震的影响，得出紫坪铺水库蓄水使得汶川地震发震断层更加危险，但1.0 kPa的库仑应力变化是否可以触发大地震的发生需要进一步研究.陶玮等(2014)采用二维有限单元法模拟紫坪铺水库蓄水引起的载荷和孔隙应力变化，得到在震源深度“ΔCFS增长值超过2~25 kPa，对汶川地震的发生有较强的促进乃至触发作用”.
可以看出，国内外不少研究组开展了紫坪铺水库是否诱发汶川地震的定量分析，提升了水库地震触发机制的认识.但是，在计算结果和结论上存在很大的分歧，有的认为紫坪铺水库蓄水引起的库仑应力可达50 kPa，与汶川地震的关系不能被忽视(雷兴林等，2008; Ge et al., 2009; Lei，2010)；有的认 为水库蓄水引起的库仑应力仅为-1.0 kPa 到10.0 kPa，与汶川地震没有物理联系(Deng et al., 2010; Gahalaut et al., 2010).那么，如何看待这些研究结果？导致计算结果不同的因素有哪些？引起的误差有多大？
2010年3月，Kerr和Stone(2010)再次在Science上发表报道“Two Years Later，New Rumblings Over Origins of Sichuan Quake”，这一争论仍在继续.目前地震学界公认的大于6.0级的水库地震全球迄今只有4个，包括我国1962年的新丰江水库地震，而这次争议的汶川地震震级高达MS8.0.无怪这一争议会成为科技界的一个热点.中国科学工作者，理应对这样一个发生在中国土地上而又引起全球关注的尖锐的问题，向世界做出一个具有说服力的分析、提出在目前科学水平上能够给出的最好答案.本文拟从水库蓄水数值计算原理出发，分析目前紫坪铺水库蓄水对汶川地震影响库仑应力变化计算文献中出现纷争的原因，找出孔隙弹性定量计算中的关键影响因素，通过数值试验了解资料其不确定性的影响.
2 水库蓄水定量分析基本原理越来越多的库区地震资料数据表明水库地震的发生与水库的蓄放水有密切关系，且不同的水库运行过程中表现出不同的响应特征(Gough and Gough.，1970; Simpson，1976; Simpson et al., 1988).一些水库初始蓄水后库区地震活动迅速增加，水库地震的因果关系非常直接明显，且具有很好的相关性，例如，我国新丰江水库1959年蓄水后不久，地震观察台便开始记录到库区的地震活动(丁原章，1989)；广西龙滩水库蓄水后地震活动不断，已发3个4级以上地震(Zhou et al., 2012)；一些水库则在蓄水后几年甚至十几年地震活动明显，表现延迟响应，相关性很小，例如，Aswan水库蓄水近17年之后，距 Aswan大坝约50 km发生了ML5.7级地震(Gahalaut et al., 2012). 有些水库蓄水后则表现出这两种机制的混合响应，例如，印度Koyna水库，在蓄水后很快出现小震，经过一系列的循环载荷变化后出现M6.5的大震，且地震活动至今仍不断.在物理力学机理上，动态响应特征的力学机制主要存在三个方面：第一，水体重力直接造成的弹性载荷改变了库区地壳受力状态，改变了给定倾角和走向的断层面上的剪应力和正应力分布，很有可能增加了地震活动可能性，当然也存在着让断层活动危险性降低的可能性；第二，蓄水造成瞬间的不排水效应，岩体中孔隙间水体来不及流出，孔隙压力立即增大，这部分会直接造成地震活动性增加；第三，孔隙压力扩散效应，对于饱和岩体，水库蓄水后在库区会形成一定的水头压力，孔隙压力逐渐向四周扩散(Simpson et al., 1988; Talwani，1997; Chen et al., 2007; Talwani et al., 2007)，直至震中位置，与不排水效应一起，影响库区应力的变化.
由于人们对于地下岩石受到的应力场缺乏实测了解，不能确定哪里应力已经临近破裂，Jaeger等(2007)据库仑准则提出了库仑破裂应力的计算公式为
设水库蓄水前构造应力场σ0<sub>ij，蓄水后应力场变化了Δσij，由震源机制解可得到发震断层几何参数(走向、倾角和滑动角)，则可求出断层面法向量(n1，n2，n3)及水库蓄水引起正应力和剪应力变化量：
水库蓄水时，水体和库体重量对库区基岩而言是一种巨大的弹性荷载，使得库基岩体弹性应力场有明显地变化.这种快速响应可看成半平面/空间边界上受分布力问题.假设水库库区岩体满足弹性力学6个基本假设条件：连续密实的物体；应力、应变和位移等是连续的；均质的；各向同性，小变形假设.当岩体颗粒之间的距离远比研究区的尺寸小时，可以将水库坝体和水库库水水体视为Boussinesq问题.水库蓄水引起库区岩体弹性应力场变化可用公 式(7)(吴家龙 2001)求解弹性应力场变化的表达式.
水库蓄水时，库区岩体外界静压力变化，若孔隙水的质量不变，与外界不发生流入/出时，孔隙水的 孔隙压力瞬间变化称为不排水效应(陈颙等，2009; Wang and Manga， 2010). 对于短时间，该不排水孔隙水压可表示为
当不考虑水体的黏性系数变化、渗透率变化等因素时，水库蓄水后孔隙压力扩散过程则可表达为(Bell and Nur， 1978; Roeloffs，1988)：
水库蓄水产生的弹性载荷、瞬间不排水孔压和孔压的扩散是一个动态耦合作用过程，据Rice和Cleary(1976)、Bell和Nur(1978)、Talwani等(1984; 1997)和Roeloffs(1988)等研究，应用孔隙弹性耦合理论对水库蓄水引起库区应力场和孔隙压力场进行分析，一方面可以模拟水库水位变化产生的弹性载荷效应(瞬间弹性、不排水响应)，另一方可以模拟孔隙压力扩散过程(排水响应)，而解耦解和非耦合解可能高估或低估了孔隙压力(Roeloffs，1988)，需采用完全孔隙弹性耦合模型计算.完全孔隙弹性耦合理论满足Biot固结理论的6个基本假设条件，其控制方程基本变量为6个应力变量、6个应变变量和孔隙压力，其本构方程：
；K是饱和岩体排水体 积模量，G是剪切模量，δij是Krnecker符号，σkk=3Kεkk；ν0是流体的容积率，Kf流体体积模量，Ks固体颗粒的体积模量.对应的渗流扩散方程：
随着模拟计算技术的提升，研究学者们据对水库库区的地质背景(地质构造、水文地质、断层、岩性等)了解与分析的基础上进行二维或者三维定量模拟分析计算.然而，库区的这些地质参数却很难测定和获取，因此，在并没有充分了解的情况下，在水库地震定量计算模型建立时通常采用对这些参数进行定性假设或者是据某个剖面(钻孔)资料估测.同时，计算模型的不同对计算结果有较大的影响，甚至可能得出性质相反的结论.对汶川这类逆断层地震，由于水库产生的压力对地震的发生可能是抑制作用，而孔隙压引的是激发作用，因而，计算模型中的断层参数的选择有着决定性的作用，如在倾角较小的情况下，计算的水库蓄水产生的库仑应力是可能抑制地震发生的(Zhou and Deng， 2011).表 1给出了不同研究学者在定量分析紫坪铺水库蓄水对汶川地震的影响时采用的模型(大小、维数)、计算方法(数值/解析)、计算震源点位置、发震断层参数、孔隙压力扩散模型、介质参数以及摩擦系数和得出的库仑应力变化值.以下，对这些影响因素做一一分析.
 | 表 1 不同学者对紫坪埔水库蓄水对汶川地震震源研究模型及研究方法Table 1 Different research models and methods by scholars on the studying of the impoundment of Zipingpu reservoir to the Wenchuan earthquake |
在研究水库地震蓄水后库区应力场变化时，有的研究者建立二维模型，有的建立三维模型.由表 1可以看出，汶川地震后，雷兴林等(2008)、Deng等(2010)分别采用三维Boussinesq解析解和二维数值解来求解紫坪铺水库蓄水引起的弹性应力和孔压扩散.Ge等(2009)、周斌等(2010)则采用二维数值模拟分别求出弹性应力和孔隙压力变化.Ghaulaut等(2012)则采用三维解析解来求解.图 1给出了三维模型和二维模型的区别，可以看出，二维模型实质上是平面应变近似，将无限的线荷载代替三维有限点载荷，势必扩大了计算结果，无形中将水库的某一剖面上的水载荷和水压放大为无限长区域.针对计 算模型维数对计算结果的影响差别，课题组郑亮等(2013)将紫坪铺水库蓄水在汶川地震震源处产生的弹性力和孔隙压力大小进行比较，发现三维模型计算的应力大小仅为二维计算的1/3~1/4.因此，在精确讨论水库蓄水触发地震定量分析时需建立三维模型.
 | 图 1 二维模型与三维模型区别示意图Fig. 1 Sketch map of the difference between 3-D and 2-D |
由水库蓄水定量计算的基本原理可看出，通过解析解和数值模拟计算均可以得出水库蓄水对库区、发震断层和震源上的库仑应力变化值.解析解方法基于一些基本假设条件，经过一系列严格的数学公式推导，得出一套通用的公式，求解过程方便简洁，可以很快得出计算结果且计算数值确定.但是，采用解析解方法计算时其假设条件比较苛刻，例如，半无限空间、各向同性的均匀介质等，然而，这些假设同实际情况还是有一定的差距.数值计算方法则很好地弥补这一缺点，可以灵活地处理地形、介质不均匀性、断层、GPS、MT测量数据等问题，较方便地通过调整初始条件和边界条件等使得计算结果更加符合实际.但是，如果计算网格不够密、数值程序存在缺陷，可能造成计算结果不对.因此，在进行数值计算时，可采用解析解来标定数值计算程序和数值解.例如，在求解水库蓄水库区应力场定量分析中，孔隙水压的扩散是一个非常复杂的问题，涉及到地下构造和渗透率分布.若假定断层两壁岩石绝对不渗透(2-D模型)，应用解析解可以给出一个孔隙压力扩散的上限值：即震源和水库之间由一条直接相连通且只有这一条断层、岩石的渗透率为0，只有这一条断层可以渗透，这时就可以用解析解求出孔压扩散的上限值.
图 2(a，b)给出了采用解析方法.紫坪铺水库自蓄水至汶川地震发生，不同扩散系数下孔隙压力随着时间-深度的扩散上限值可以看出，扩散系数的选取对计算结果有较大的影响.当选取扩散系数为10.0 m2·s-1，库心至震源距为18~35 km时，孔隙压力的上限值为0.39~0.66 MPa；扩散系数为0.5 m2·s-1，上限值为1.0 kPa~0.05 MPa.然而，紫坪铺水库附近存在多条大小断层且不能认为库区周围岩石渗透率为0.水库水载荷和孔隙压力是多方向渗透，并不是只沿一条断层向震源渗透，因此，孔隙水压增高不会如上解析解计算的那么快和大.但是，给出了一个上限值，数值解大于它一定错，小于它则可能对.用该解析解可以讨论不同作者条件下的上限，判别不同作者是否合理；也可以对震源位置的影响有所了解.另外，解析解目前没有办法提供瞬间孔隙压力变化随着时间扩散，这部分只有数值解计算才能提供，但其影响较小，见图 2c，紫坪铺 水库水深100 m、扩散系数为0.2 m2·s-1时，产生的非排水孔压瞬间在浅层岩体增加较大(约2.5 kPa)，但此部分孔压随着时间向深部岩体扩散较快，在汶川地震发生时，距水库库心10 km处仅约100 Pa.
 | 图 2 不同扩散系数，紫坪铺水库自蓄水至汶川地震发生时孔隙压力随着时间-深度的上限值 (a)不同扩散系数下不同深度孔隙压力在蓄水32个月和平均100 m水深的孔隙压力扩散值;(b)库心至震源距分别为18 km和35 km时，不同孔隙压力扩散上限值;(c)瞬间不排水孔压随着时间-深度扩散(紫坪铺水库于2005年9月开始蓄水，距汶川地震发生32个月. 水库水位由752 m达到877 m，平均库水深100 m).Fig. 2 The calculated limited diffusive pore pressure of the hypocenter variations with depth and time for different diffusivity coefficient values (a)The calculated diffusive pore pressure of the hypocenter variations with depth when the water level of Zipingpu reservoir is 100 m and the storage time is 32 months;(b)The calculated diffusive pore pressure of hypocenter at the distance of 18 km and 35 km when the water level of Zipingpu reservoir is 100 m;(c)The instantaneous undrain pore pressure with depth and time(The first filling of Zipingpu reservoir was in September 2005，befor 32 months of Wenchuan earthquake occurred. The water levels from 752 m to 877 m and average depth is 100 m). |
解析方法局限于简单条件，一旦考虑复杂的地质背景将变得非常复杂，甚至无解.因此，在复杂情况下，须采用数值计算方法，将更多的地质资料、地质背景加以考虑.紫坪铺水库位于青藏高原与成都平原交界处东部，靠近龙门山断裂带，处于一个复杂的地质构造体系中.在综合探讨紫坪铺水库蓄水对汶川地震的影响，须采用数值计算方法将更多的地质资料，例如，岩体介质分层、断层性质、库区地形高程、水位动态变化等加以考虑，但需要保证一定的网格精度，减小数值误差.鉴于已有研究组采用的数值计算模型(Ge et al., 2009; Deng et al., 2010)，我们给出了2-D模型中孔隙水压扩散数值计算中不同网格密度计算结果与精确的解析解结果对比，见图 3.可以看出，网格粗细对数值计算结果有一定的影响，特别是在水平剖面上.在85 km×40 km的2-D模型中，当网格小于1.0 km ×1.0 km时，18 km深度数值绝对误差最大为0.15 kPa，是可以接受的，而网格精度差时引起的数值误差将不能忽视.因此，在详细分析谈论水库蓄水引起库区应力场变化时，采 用数值计算须满足一定的网格精度以减少计算误差.
 | 图 3 采用数值计算时采用不同网格密度对计算结果的影响 (a)网格精度为1 km × 1 km，紫坪铺水库水深100 m、蓄水32个月时孔隙压力扩散;(b)采用不同网格精度(2.5 km×2.5 km，1.0 km ×1.0 km，0.5 km×0.5 km，0.25 km×0.25 km)时，横剖面H-H′(深度18 km处)孔压扩散与解析解对比;(c)采用不同网格精度(2.5 km×2.5 km，1.0 km×1.0 km，0.5 km×0.5 km，0.25 km×0.25 km)时，纵剖面V-V′(库心正下方)孔压扩散与解析解对比.Fig. 3 Effects of the grids size on the calculation results by numerical method (a)Distribution of diffusive pore pressure and the grid size is 1 km × 1 km，when the water level of Zipingpu reservoir is 100 m and the storage time is 32 months，the width of the reservoir is 5 km.(b)The diffusive pore pressure of profile H-H′ variation with distance for different grid sizes(2.5 km×2.5 km，1.0 km×1.0 km，0.5 km×0.5 km，0.25 km×0.25 km).(c)The diffusive pore pressure of profile V-V′ variation with distance for different grid sizes(2.5 km×2.5 km，1.0 km×1.0 km，0.5 km×0.5 km，0.25 km×0.25 km). |
汶川地震后，各个研究学者计算了紫坪铺水库蓄水后震源处的库仑应力变化(ΔCFS)，进而探讨水库对大地震的影响.然而，据公式(4—6)，在计算中，发震断层参数、震源和库心位置的选取对计算结果有着较大的影响.从表 1中，可以看出不同研究者在计算ΔCFS时，震源位置、采用用断层参数各有所不同.椐Boussinesq′s解，弹性应力场随着深度呈现出2次方衰减.由图 2a可以看出，孔隙压力扩散与震源和库心之间的距离有着密切关系.汶川地震后，美国地质调查局(USGS－NEIC)、哈佛大学CMT以及国内台网中心张勇等(2008)、吕坚等(2008)、王卫民等(2008)分别给出汶川地震震源机制解.然而，由于各个研究机构采用的计算方法不同得出的结果不尽相同，特别是震源深度，见表 2和图 4，而 在计算中各个研究组选取紫坪铺水库库心也不一致.为此，我们将估算断层参数和震源机制对计算结果的影响.
| 表 2 不同机构给出的汶川地震发震断层参数、震源位置及其椐水库库心距离Table 2 Fault planes derived from different institutes and the distance between the Zipingpu Reservoir and the hypocenter of Wenchuan earthquake |
 | 图 4 不同机构给出的汶川地震震源位置和水库库心点A、B和C分别为在计算点载荷时紫坪铺水库库心的选取：点A是紫坪铺大坝的地理位置；点B是紫坪铺水库几何中心；点C为整 个水库的几何中心.对应的坐标值分别是(103.5，31.016)，(103.54，31.025)和(103.574，31.04).震源球a、b、c和d分别表示不同机构给出汶川地震震源位置：震源球a是美国地质调查局(USGS-NEIC)给出；震源球b是吕坚等(2008)给出； 震源球c是台网中心张勇等(2008)给出和震源球d是王卫民等(2008)给出.Fig. 4 Distribution of the epicenters and fault zones in the Zipingpu reservoir area The A、B and C are the location of the point sources. Point A is nearly the location of the Zipingpu dam. Point B is closed to the center of the reservoir area. And point C is the center of whole reservoir area. The corresponding locations are(103.5，31.016)，(103.54，31.025)，(103.574，31.04). The focal spheres(a、b、c and d)are from different research institutes. The focal sphere a is given from USGS-NEIC; the focal sphere b is from lv，et al.(2008); the focal sphere c is from Zhang，et al.(2008); the focal sphere a is from Wang，et al.(2008) |
紫坪铺水库蓄水1.12 km3，正常水位877 m(http://en.wikipedia.org/ wiki/Zipingpu_Dam)，基于Boussinesq解，分别选取点A，B和C作为点载荷点，我们给出3-D情况下不同震源参数下紫坪铺水库蓄水后库区弹性应力场变化，见表 3，可以看出计算得出的弹性应力场产生的库仑应力变化有着较大差别.除CMT提供的矩心矩张量震源机制外，选取不同点载荷点和震源机制情况下，紫坪铺水库 蓄水产生的弹性库仑应力变化为-0.3～-4.9 kPa，相差达10倍以上.汶川地震发生前，紫坪铺水库水位减小至821 m，我们同样计算了此部分放水产生的弹性应力变化产生的库仑应力变化，见表 4，放水使得震源处的库仑应力增加，选取B点作为库心时，不同震源参数下ΔCFS范围为0.7～2.2 kPa，相差达到3倍以上.进一步计算了不同震源参数下孔隙压力扩散，叠加上弹性应力场变化可得出总库仑应力变化，见表 5，不同震源参数下的ΔCFS相差达到2~7倍以上.
| 表 3 紫坪铺水库蓄水不同点载荷下不同震源机制弹性应力场变化Table 3 Elastic stresses changes induced by the impoundment of Zipingpu reservoir |
| 表 4 坪铺水库水位由877 m下降至821 m时，不同震源机制弹性应力场变化(单位：Pa)Table 4 Elastic stresses changes due to the discharge of Zipingpu reservoir(water levels from 877 to 821 m; Point B as the reservoir gravity center; the friction coefficient is 0.6)(Units: Pa) |
| 表 5 紫坪铺水库水位由877 m下降至821 m时，不同震源参数下库仑应力变化Table 5 ΔCFStotal induced by the impoundment of Zipingpu reservoir |
水库蓄水后，库水位增高，孔隙压力增大且以水库库体形态进行体源扩散.虽然2-D模型计算结果同3-D模型相比存在一定的误差，夸大计算结果.但为了更快、更方便地分析库区应力和孔压扩散，常常会建立2-D模型.然而，在建立2-D模型时，不同研究组则可能选取不同的剖面进行孔隙压力沿断层面扩散的模拟计算.例如，汶川大地震后，Ge等(2009)选取了接近水库库尾且垂直于紫坪铺水库走向和映秀—北川断裂的垂直剖面；周斌等(2010)则在接近水库中心位置、垂直于水库走向的垂直剖面；而Deng等(2010)选取的是沿着紫坪铺水库走向且通过震源点的斜面，见图 5.显然，在选取的不同剖面上水库紫坪铺水库蓄水产生的应力场、孔压场变化不同，且距离震源的位置也不同.同时，在各自选取的2-D剖面计算中若依然按照震源机制提供的震源深度或者震源在该剖面上的水平投影均会产生一定的误差.此外，在数值计算模拟水库蓄水后库区应力场变化时，最理想状态是将库区的网格划分无穷小，据水位和库底高程值给出每一点的应力载荷和初始孔隙压力.在紫坪铺水库2-D模拟中，Ge等(2009)和Deng等(2010)在各自2-D模型上均假定紫坪铺水库宽5 km，在这5 km库区前者给定应力载荷和初始孔隙压力由水库蓄水100 m水深产生，而后者给定为50 m水深(5 bar)，可以看出即使两者选取的网格精度、扩散系数、扩散模型和震源位置一样的条件下，两者计算得出的应力值也会有较大的区别，前者得出的计算结果要大于后者.因此，在详细分析水库蓄水对库区应力场影响时需进行3-D体扩散模拟.
 | 图 5 不同研究组给出的紫坪铺水库2-D模型Fig. 5 The 2-D profiles of in different models for the Zipingpu reservoir |
在孔隙压力计算中孔隙压力扩散系数是最重要的参数，由图 2也可以看出在不同扩散系数孔隙压力在同一深度处相差几百倍.目前，对于水库地震扩散系数(渗透率)的选取众多学者有不同的建议.一方面是据已有的地震进行统计分析，例如，Talwani和Acree(1984)据水库蓄水后地震活动面积与地震活动滞后时间关系得出了地震扩散系数，其大小在1～10 m2·s-1 之间.另一方面，通过室内实验室测定岩体的渗透率，再推导孔压扩散系数，例如，龚钢延和谢原定(1989; 1991)通过实验研究得出完整花岗岩渗透率1.8×10-18m2，破裂花岗岩的渗透率1.8×10-15m2，对应的孔隙扩散系数仅为0.2 m2·s-1，两者相差还是比较大的.在紫坪铺水库计算上，国内 外研究组在计算中采用的孔隙压力扩散系数也不相同，见表 1，取值范围相当广泛，0.001~100 m2·s-1.因此，对扩散系数的取舍，需要结合现有资料基础，通过野外工作和室内试验综合分析取得更精确地认识.
4 结论与讨论 4.1 为什么不同研究组给出不同的计算结果？汶川地震后，众多学者对距震中近10 km的紫坪铺水库是否可以触发/诱发MS8.0地震发生进行激烈的讨论.本文分析总结了国内外关于紫坪铺水库对汶川地震影响研究的基本情况，针对各个研究者的研究方法、计算模型和模拟方法及控制性参数的选取等方面进行分析，系统地分析和评述得出紫坪铺水库对汶川地震的主要影响因素，以期解决目前水库蓄水定量计算中出现的争议.初步得出：①有些属于方法的缺陷：例如，采用二维模拟而不采用三维模拟，会得到三倍以上夸大的结果；仅考虑断层渗透率(把岩体渗透率视为无穷大)或仅考虑均匀各向同性的岩体渗透率(忽视断层渗透率)，均具有片面性；弹性载荷解析解无法考虑弹性载荷作用下瞬间孔隙弹性介质中产生的孔隙压力及其随后的扩散对孔隙压力的贡献.②有些属于模型参数取舍的不同：例如震源位置和深度，深度相差一倍，孔隙压力传递所需时间会增加4倍；震源机制解断层走向倾角的差异，会显著影响库仑应力大小计算结果，可到达2~7倍；蓄/放水加压量值和空间分布处理的不同；岩体弹性性质和孔隙弹性性质取值的不同；岩体和断层扩散率取值的不同等；这些对结果都会有程度不同的影响.③有些属于数值计算不同，例如，计算程序网格精度及计算精度等问题的处理上存在不严格之处.④水库蓄水产生的弹性载荷、瞬间不排水孔压和孔压的扩散是一个动态耦合作用过程，有些是简单地将2-D孔隙压力数值解同3-D弹性应力场解析解简单叠加，或采用孔隙弹性耦合度解耦解和非 耦合解可能高估或低估了孔隙压力(Roeloffs，1988).
4.2 裂隙岩体扩散率的进一步研究必要性经历了构造运动，岩体中存在着大量的断层、节理、裂隙等规模不同的结构面.这些软弱的结构面不仅仅控制着岩体力学特征，同时，也是岩体中地下水等流体的通道，是孔隙压力快速传递扩散的通道.在水库蓄水引起库区孔隙压力变化定量化计算中，岩体渗透率/孔隙压力扩散系数是最重要的一个参数.在一定程度上扩散系数的选取多少决定了计算结果大小.目前，对于扩散率的选取一方面是对已经发生确定的水库地震进行的定性统计，另一方面在实验室测定的往往是完整岩石的渗透率/孔压扩散系数，两者均存在一定的局限性.因此，对岩体渗透率/孔隙压力扩散系数研究是水库地震研究中一个必不可少的部分.
我们知道岩体中的节理裂隙系统往往由几类产状不同的节理组合而成，而同组的裂隙在形成原因、分布特点和工程性质方面都具有相似性.因此，可以在野外测得节理的倾向和倾角，虽然变化有时候很大，但可以进一步分组分类，进行一些分布函数的拟合，或者采用蒙特卡洛等方法建立一概率模型用伪随机数作近似抽样进行统计整理得出，见图 6，从而可以采用数值模拟方法来计算水库库区孔隙压力扩散系数.另一方面，在野外采集库区岩样，应用X-ray层析与LBM方法结合，测定非完整岩石(松软、破碎)等难以取得大试样的渗透率.
 | 图 6(a)采用数值模拟方法计算渗透率/孔隙扩散系数;(b)、(c)来自参考文献(Jing等，2010; 孙玉军等，2013)Fig. 6 The sketch map of using numerical simulation method to calculate permeability or diffusion coefficient;(b)，(c)from Jing et al., 2010 and Sun et al., 2013) |
水库蓄水后水体对库底岩体应力状态的改变是一个非常复杂的过程，涉及到库区的背景构造应力场、岩性、构造的分布及其参数等问题，往往这些地质情况决定了水库地震何时何地产生.例如，马文涛等(2012)使用灰色聚类方法，对已发生水库地震库区的岩性、地应力状态等进行统计，对于可溶性岩、挤压/拉张区域更易发生震级较大的水库地震.广西龙滩水库自2006年蓄水现已发生3个4级以上地震，而水库地震的发生恰恰处于电阻降低的区域(詹艳等，2012).Zhou等(2012)通过对龙滩水库库区的VP/VS的研究发现库区地震活动主要集中在低VP区域.因此，对水库地震触发机制的研究中需要综合考虑水库库区实际地球物理观测数据、地质条件、构造条件、蓄放水历史记录、GPS数据、MT数据和库区地震目录约束条件，而有限元数值模拟方法可以将这些约束条件灵活地引入，并得到高分辨率的数值结果.简单的一维或二维模型，通过解析解的方式已经无法满足实际研究的需要了.此外，若能够详尽考虑实际观测结果，不采用高分辨率的有限元数值建模也是不可能的.在这种情况下，部分区域的分辨率需要达到几十米到几米的分辨率，而绝大部分区域需要公里分辨率的网格.综合下来，一个水平和深度方向都在几十公里范围的库区建模，至少需要100~300万有限元网格才能达到精度需求.因此，在进行定量分析时需要将这些地质数据加入到地质模型中，需要通过高分辨率数值模拟计算来实现，从而揭示更为合理的库区弹性应力、孔隙压力和应变能的动态变化规律.
4.4 是否有库仑应力变化触发阈值？关于多大量级的库仑应力变化可以触发地震已有不少研究.例如，Stein等(1983; 1992)计算主震后库仑应力变化，得出库仑应力变化对余震的影响的典型值在0.1~1.0 MPa.Toda等(1998)研究Mw6.9日本西南部Kobe主震应力触发余震空间分布的关系，发现余震集中在库仑应力变化0.02~0.1 MPa区域.Bell和Nur(1978)通过简化的二维模型估算了200 m深的奥洛维尔水库对正断裂引起的库仑应力变化量级为0.08~0.14 MPa，且孔隙压力超过0.35 MPa地震活动明显增加.沈立英等(1992)应用二维孔隙弹性解耦模型得出新丰江水库MS6.1级地震由于新丰江水库蓄水引起震源处库仑应力变化为0.3 MPa.刘桂萍和傅征详(2009)计算了我国海原地震对古浪地震、唐山地震区域库仑应力变化，得出库仑应力变化大小约为0.01 MPa的数量级，可触发区域地震的发生.Pauchet等(1999)对法国ML5.2地震之后的余震主要集中在库仑应力变化在0.2 MPa区域.King等(1994)认为在一定条件下0.01 MPa量级微小力变化可明显影响地震活动.之后，研究学者对各个大震的震后应力变化，验证了正的库仑应力变化能使断层上的地震活动性增加，而在库仑应力变化为负的区域地震活动性降低(Hardebeck et al., 1998)，且应力触发阈值在0.01 MPa量级(Harris，1998).
更小量级的库仑应力变化是否可以触发地震？Ziv和Rubin(2000)发文探讨对于地震应力触发问题，统计California地区63个M>4.5地震，主要得 出对于该区库仑应力变化量级小于0.01 MPa，在选取不同的摩擦系数后库仑应力变化量级小于0.01 MPa 及0.001 MPa的地震数目仍然占到了65%.另外，不少研究学者指出潮汐应力触发地震的库仑应力变化虽然量级(0.001 MPa)很小，但也不排除其触发地震的可能性(Heaton，1975; Vidale et al., 1998; Perfettini and Schmittbuhl， 2001).
“紫坪铺水库蓄水是否能够触发汶川大地震的发生？”，鉴于目前的研究成果，库仑应力变化在kPa量级，尚不能排除触发的可能性，但得出的蓄水震源处的库仑应力变化太低，在背景构造应力场不明确的情况下，也不能确定一定有联系.从断裂力学破裂准则角度分析，在应力临界点处只要很小的应力波动即可触发断裂.那么水库地震作为一特殊类型的地震，到底可估算的多大量级的库仑应力变化可以触发水库地震？程惠红等(2012)对我国公认的1962年新丰江MS6.1级水库地震造成的库仑应力变化进行量化估算，得出新丰江水库蓄水引起库区和其周围断层上库仑应力变化量级为0.001 MPa.鉴于水库地震的复杂性，仅仅对少数一、两个水库地震的研究模拟不能真正认识该问题的全貌.因此，需进一步选取不同类型的中国和世界上已发生过地震的水库作为研究对象，进行定量分析水库蓄水后引起库区应力场、应变和库仑应力变化，计算分析在水库运行过程中不同类型的库区应力场特征.同时，水库地震的发生是在水与断层、岩体相互作用的动态演化过程中发生的，因此，在今后工作中需要考虑水岩动态作用过程，进而可以探讨不同机制下水库地震触发机制特点，丰富水库地震触发物理机制的定量认识.
| [1] | Bell M L, Nur A. 1978. Strength changes due to reservoir-induced pore pressure and stresses and application to lake Oroville. J. Geophysics. Res., 83(B9):4469-4483. |
| [2] | Chen H Q, Xu Z P, Lee M. 2008. Wenchuan Earthquake and seismic safety of large dams. Shuili Xuebao (in Chinese), 39(10):1158-1167. |
| [3] | Chen L, Talwani P. 2001. Mechanism of initial seismicity following impoundment of the Monticello Reservoir, South Carolina. Bulletin of the Seismological Society of America, 91(6):1582-1594. |
| [4] | Chen L, Talwani P. 1998. Reservoir-induced seismicity in China. Pure and Applied Geophysics, 153(1):133-149. |
| [5] | Chen P, Wang J D. 2007. Wavelet estimation of the diffusion coefficient in time dependent diffusion models. Science in China Series A:Mathematics, 50(11):1597-1610. |
| [6] | Chen Y. 2009. Did the reservoir impoundment trigger the Wenchuan earthquake?. Science in China Series D:Earth Sciences, 52(4):431-433. |
| [7] | Chen Y, Huang T F, Liu N R. 2009. Rock Physics (in Chinese). Hefei:Press of University of Science and Technology of China. |
| [8] | Cheng H H, Zhang H, Zhu B J, et al. 2012. Finite element investigation of the poroelastic effect on the Xinfengjiang Reservoir-triggered earthquake. Science in China Series D:Earth Sciences, 55(12):1942-1952. |
| [9] | Cheng W Z, Zhang Z W, Ruan X. 2010. Analysis of seismicity and its causes in Zipingpu reservior region at each water storage stage. Progress in Geophys. (in Chinese), 25(3):759-767, doi:10.3969/j.issn.1004-2903.2010.03.003. |
| [10] | Deng K, Zhou S, Wang R, et al. 2010. Evidence that the 2008 Mw7.9 Wenchuan Earthquake Could Not Have Been Induced by the Zipingpu Reservoir. Bull. Seism. Soc. Am., 100(5B):2805-2814. |
| [11] | Ding Y Z. 1989. The Reservoir Induced Earthquake (in Chinese). Beijing:Seismological Press. |
| [12] | Gahalaut K, Gahalaut V K. 2010. Effect of the Zipingpu reservoir impoundment on the occurrence of the 2008 Wenchuan earthquake and local seismicity. Geophysical Journal International, 183(1):277-285. |
| [13] | Gahalaut K, Hassoup A. 2012. Role of fluids in the earthquake occurrence around Aswan reservoir, Egypt. Journal of Geophysical Research, 117(B2):B02303, doi:10.1029/2011JB008796. |
| [14] | Ge S M, Liu M, Lu N, et al. 2009. Did the Zipingpu Reservoir trigger the 2008 Wenchuan earthquake?. Geophysical Research Letters, 36(20):L20315, doi:10.1029/2009GL040349. |
| [15] | Gong G Y, Xie Y D. 1989. A review of laboratory investigation of rock permeability. Chinese Journal Rock Mechanics and Engineering (in Chinese), 8(3):219-227. |
| [16] | Gong G Y, Xie Y D. 1991. Research on the diffusion of pore fluid and in-situ hydraulic diffusivity in the epicentral region of Xinfengjiang reservoir earthquakes. Acta Seismologica Sinica (in Chinese), 13(3):364-371. |
| [17] | Gough D I, Gough W L. 1970. Stress and deflection in the lithosphere near Lake Kariba—I. Geophys. J. Roy., 21(1):65-78. |
| [18] | Gupta H K. 2002. A review of recent studies of triggered earthquakes by artificial water reservoirs with special emphasis on earthquakes in Koyna, India. Earth-Science Reviews, 58(3-4):279-310. |
| [19] | Hardebeck J L, Nazareth J J, Hauksson E. 1998. The static stress change triggering model:Constraints from two southern California aftershock sequences. J. Geophys. Res., 103(B10):24427-24437. |
| [20] | Harris R A. 1998. Introduction to special section:Stress triggers, stress shadows, and implications for seismic hazard. J. Geophys. Res., 103(B10):24347-24358. |
| [21] | Heaton T H. 1975. Tidal triggering of earthquakes. Geophys. J. Int., 43(2):307-326. |
| [22] | Jaeger J C, Cook N G W, Zimmerman R. 2007. Fundamentals of Rock Mechanics. New York:John Wiley & Sons. |
| [23] | Jing H, Yu X, Zhang H, et al. 2010. Numerical analysis on thermal conductivity of poly-mineral rock. Earthquake Science, 23(3):223-232. |
| [24] | King G C P, Stein R S, Lin J. 1994. Static stress changes and the triggering of earthquakes. Bulletin of the Seismological Society of America, 84(3):935-953. |
| [25] | Klose C D. 2012. Evidence for anthropogenic surface loading as trigger mechanism of the 2008 Wenchuan earthquake. Environmental Earth Sciences, 66(5):1439-1447. |
| [26] | Lei X L. 2010. Possible roles of the Zipingpu Reservoir in triggering the 2008 Wenchuan earthquake. Journal of Asian Earth Sciences, 40(4):844-854. |
| [27] | Lei X L, Ma S L, Wen X Z, et al. 2008. Integrated analysis of stress and regional seismicity by surface loading—A case study of Zipingpu reservoir. Seismology and Geology (in Chinese), 30(4):1046-1064. |
| [28] | Liu G P, Fu Z X. 2001. A study on the great gulang earthquake triggered probably by static stress change resulting from the great haiyuan earthquake. Chinese J. Geophys. (in Chinese), 44(z1):107-115, doi:10.3321/j.issn:0001-5733.2001.z1.014. |
| [29] | Lü J, Su J R, Jin Y K, et al. 2008. Discussion on relocation and seismo-tectonics of the Ms8.0 Wenchuan earthquake sequences. Seismology and Geology, 30(4):917-925. |
| [30] | Ma W T, Xu X W, Yu G H, et al. 2012. Assess the reservoir- induced seismic hazard in the hubei section of the three gorges reservoir using gray clustering method. Seismology and Geology (in Chinese), 34(4):726-738. |
| [31] | Ma W T, Xu C P, Zhang X D, et al. 2011. Study on the relationship between the reservoir- induced seismicity at zipingpu reservoir and the Ms8.0 wenchuan earthquake. Seismology and Geology (in Chinese), 33(1):175-190. |
| [32] | Pauchet H, Rigo A, Rivera L, et al. 1999. A detailed analysis of the February 1996 aftershock sequence in the eastern Pyrenees, France. Geophys. J. Int., 137(1):107-127. |
| [33] | Perfettini H, Schmittbuhl J. 2001. Periodic loading on a creeping fault:Implications for tides. Geophysical Research Letters, 28(3):435-438. |
| [34] | Rice J R, Cleary M P. 1976. Some basic stress diffusion solutions for fluid-saturated elastic porous media with compressible constituents. Reviews of Geophysics, 14(2):227-241. |
| [35] | Roeloffs E A. 1988. Fault stability changes induced beneath a reservoir with cyclic variations in water level. J. Geophys. Res., 93(B3):2107-2124. |
| [36] | Shen L Y, Chang B Q, Wu M B. 1992. Application of the stress-pore pressure coupled theory of porous medium to the research of inducing earthquakes at Xinfengjiang Reservoir. Earthquake Research China (in Chinese), 8(3):36-43. |
| [37] | Simpson D W. 1976. Seismicity changes associated with reservoir loading. Eng. Geol., 10(2-4):123-150. |
| [38] | Simpson D W, Leith W S, Scholz C H. 1988. Two types of reservoir-induced seismicity. Bull. Seism. Soc. Am., 78(6):2025-2040. |
| [39] | Stein R S, Lisowski M. 1983. The 1979 Homestead Valley earthquake sequence, California:Control of aftershocks and postseismic deformation. J. Geophys. Res., 88(B8):6477-6490. |
| [40] | Stein R S, King G C P, Lin J. 1992. Change in failure stress on the southern San Andreas fault system caused by the 1992 magnitude=7.4 Landers earthquake. Science, 258(5086):1328-1332. |
| [41] | Sun Y J, Kim I K, Zhang H, Shi Y L. 2013. Simulation of tensile test of the steel fiber reinforced concrete with the finite element method. Earthquake(in Chinese), 33(4):145-152. |
| [42] | Sun Y J, Zhang H, Dong S W, et al. 2012. Study on effect of the Zipingpu reservoir on the occurrence of the 2008 Wenchuan earthquake based on a 3D-poroelastic model. Chinese J. Geophys. (in Chinese), 55(7):2353-2361, doi:10.6038/j.issn.0001-5733.2012.07.020. |
| [43] | Talwani P. 1997. On the nature of reservoir-induced seismicity. Pure and Applied Geophysics, 150(3-4):473-492. |
| [44] | Talwani P, Acree S. 1984. Pore pressure diffusion and the mechanism of reservoir-induced seismicity. Pure and Applied Geophysics, 122(6):947-965. |
| [45] | Talwani P, Chen L, Gahalaut K. 2007. Seismogenic permeability, ks. J. Geophys. Res., 112(B7):B07309, doi:10.1029/2006JB004665. |
| [46] | Tao W, Masterlark T, Shen Z K, et al. 2014. Triggering effect of the Zipingpu Reservoir on the 2008 Mw7.9 Wenchuan, China, Earthquake due to poroelastic coupling. Chinese J. Geophys. (in Chinese), 57(10):3318-3331, doi:10.6038/cjg20141019. |
| [47] | Toda S, Stein R S, Reasenberg P A, et al. 1998. Stress transferred by the1995 Mw=6.9 Kobe, Japan, shock:effect on aftershocks and future earthquake probabilities. J. Geophys. Res., 103(10):24543-24565. |
| [48] | Vidale J E, Agnew D C, Johnston M J S, et al. 1998. Absence of earthquake correlation with Earth tides:An indication of high preseismic fault stress rate. J. Geophys. Res., 103(B10):24567-24572. |
| [49] | Wang C Y, Manga M. 2010. Earthquakes and Water. New York:Springer. |
| [50] | Wang W M, Zhao L F, Li J, et al. 2008. Rupture process of the Ms8.0 Wenchuan earthquake of Sichuan, China. Chinese J. Geophys.(in Chinese), 51(5):1403-1410. |
| [51] | Wu J L. 2001. Mechanics of Elasticity (in Chinese). Beijing:China Higher Education Press. |
| [52] | Zhang B, Shi Y L. 2010. A discussion on the influences of Zipingpu reservoir on the stability of faults in the neighborhood. Journal of the Graduate School of the Chinese Academy of Sciences, 27(6):754-760. |
| [53] | Zhang Y, Feng W P, Xu L S, et al. 2009. Spatio-temporal rupture process of the 2008 great Wenchuan earthquake. Science in China Series D:Earth Sciences, 52(2):45-154. |
| [54] | Zhang Y, Wang L F, Wang J J, et al. 2012. Electromagnetic survey of the seismogenic structures beneath the Longtan reservoir in Guangxi Province. Chinese J. Geophys. (in Chinese), 55(4):1400-1410, doi:10.6038/j.issn.0001-5733.2012.04.036. |
| [55] | Zheng L, Zhang H, Sun Y J, et al. 2013. Comparison of 2D and 3D finite element simulations in the research of reservoir induced earthquakes. Earthquake (in Chinese), 33(4):162-171. |
| [56] | Zhou B, Xue S F, Deng Z H, et al. 2010. Relationship between the evolution of reservoir-induced seismicity in space-time and the process of reservoir water body load-unloading and water infiltration—A case study of Zipingpu reservoir. Chinese J. Geophys. (in Chinese), 53(11):2651-2670, doi:10.3969/j.issn.0001-5733.2010.11.013. |
| [57] | Zhou L Q, Zhao C P, Chen Z L, et al. 2012. Three-dimensional VP and VP/VS structure in the Longtan reservoir area by local earthquake tomography. Pure and Applied Geophysics, 169(1-2):123-139. |
| [58] | Zhou S, Deng K. 2011. Reply to "Comment on 'Evidence that the 2008 Mw7.9 Wenchuan Earthquake could not have been induced by the Zipingpu Reservoir' by Kai Deng, Shiyong Zhou, Rui Wang, Russell Robinson, Cuiping Zhao, and Wanzheng Cheng" by Shemin Ge. Bull. Seism. Soc. Am., 101(6):3119-3120. |
| [59] | Zhu B, Liu C, Shi Y L, et al. 2011. Apllication of flow driven pore-network crack model to Zipingpi reservoir and Longmenshan slip. Science China (Physics, Mechanics & Astronomy), 54(8):1532-1540. |
| [60] | Ziv A, Rubin A M. 2000. Static stress transfer and earthquake triggering:No lower threshold in sight?. J. Geophys. Res., 105(B6):13631-13642. |
| [61] | 陈厚群, 徐泽平, 李敏. 2008. 汶川大地震和大坝抗震安全. 水利学报, 39(10): 1158-1167. |
| [62] | 陈颙. 2009. 汶川地震是由水库蓄水引起的吗?. 中国科学D: 地球科学, 39(3): 257-259. |
| [63] | 陈颙, 黄庭芳, 刘恩儒. 2009. 岩石物理学. 合肥: 中国科学技术大学出版社. |
| [64] | 程惠红, 张怀, 朱伯靖等. 2012. 新丰江水库触发地震的孔隙弹性耦合有限元模拟. 中国科学D: 地球科学, 42(6): 905-916. |
| [65] | 程万正, 张致伟, 阮祥. 2010. 紫坪铺水库区不同蓄水阶段的地震活动及成因分析. 地球物理学进展, 25(3): 759-767, doi: 10.3969/j.issn.1004-2903.2010.03.003. |
| [66] | 丁原章. 1989. 水库诱发地震. 北京: 地震出版社. |
| [67] | 龚钢延, 谢原定. 1989. 岩石渗透率变化的实验研究. 岩石力学与工程学报, 8(3): 219-227. |
| [68] | 龚钢延, 谢原定. 1991. 新丰江水库地震区内孔隙流体扩散那与原地水力扩散率的研究. 地震学报, 13(3): 364-371. |
| [69] | 雷兴林, 马胜利, 闻学泽等. 2008. 地表水体对断层应力与地震时空分布影响的综合分析——以紫坪铺水库为例. 地震地质, 30(4): 1046-1064. |
| [70] | 刘桂萍, 傅征祥. 2001. 海原大地震对古浪大地震的静应力触发研究. 地球物理学报, 44(增刊): 107-115, doi: 10.3321/j.issn:0001-5733.2001.z1.014. |
| [71] | 吕坚,苏金蓉,靳玉科等. 2008. 汶川8.0级地震序列重新定位及其发震构造初探. 地震地质,30(4):917-925. |
| [72] | 马文涛, 徐锡伟, 于贵华等. 2012. 使用灰色聚类方法评估长江三峡水库湖北不同库段水库诱发地震的震级上限. 地震地质, 34(4): 726-738. |
| [73] | 马文涛, 徐长朋, 张新东等. 2011. 紫坪铺水库与汶川地震关系的讨论. 地震地质, 33(1): 175-190. |
| [74] | 沈立英, 常宝琦, 吴名彬. 1992. 弹性多孔介质中应力-孔隙压耦联理论在新丰江水库地震研究中的应用. 中国地震, 8(3): 36-43. |
| [75] | 孙玉军, Kim I K, 张怀等.2013. 利用有限元方法模拟力学拉伸试验——以钢纤维混凝土为例. 地震, 33(4): 145-152. |
| [76] | 孙玉军, 张怀, 董树文等. 2012. 利用三维孔隙弹性模型探讨紫坪铺水库对汶川地震的影响. 地球物理学报, 55(7): 2353-2361, doi: 10.6038/j.issn.0001-5733.2012.07.020. |
| [77] | 陶玮, Masterlark T, 沈正康等. 2014. 紫坪铺水库造成孔隙弹性耦合变化及其对2008年汶川地震触发作用. 地球物理学报, 57(10): 3318-3331, doi: 10.6038/cjg20141019. |
| [78] | 王卫民, 李连锋, 李娟等. 2008.四川汶川8.0级地震震源过程.地球物理学报.51(5):1403-1410. |
| [79] | 吴家龙. 2001. 弹性力学. 北京: 高等教育出版社. |
| [80] | 张贝, 石耀霖. 2010. 紫坪铺水库对附近断层稳定性影响的探讨. 中国科学院研究生院学报, 27(6): 754-760. |
| [81] | 张勇, 冯万鹏, 许力生等. 2008. 2008年汶川大地震的时空破裂过程. 中国科学D: 地球科学, 38(10): 1186-1194. |
| [82] | 詹艳, 王立凤, 王继军等. 2012. 广西龙滩库区深部孕震结构大地电磁探测研究. 地球物理学报, 55(4): 1400-1410, doi: 10.6038/j.issn.0001-5733.2012.04.036. |
| [83] | 郑亮, 张怀, 孙玉军等. 2013. 水库触发地震研究中二维与三维有限元模拟结果比较. 地震, 33(4): 162-171. |
| [84] | 周斌, 薛世峰, 邓志辉等. 2010. 水库诱发地震时空演化与库水加卸载及渗透过程的关系——以紫坪铺水库为例. 地球物理学报, 53(11): 2651-2670, doi: 10.3969/j.issn.0001-5733.2010.11.013. |