引言

建筑物的抗震设防需要尽可能地掌握未来可能遭遇的大地震引起的强地震动信息,但中国地震记录较为短缺,需要靠模拟地震动来弥补真实地震记录的不足。地震动模拟常用的方法有理论格林函数方法、经验格林函数方法(Irikura, 1983, 1986, 1989, 2000, 2006Irikura等,1994, 2004)和随机方法(Boore,2003)。理论格林函数模拟地震动虽然在理论上比较严密,但在实际应用过程中有2大缺点,一是理论格林函数模拟地震动必须掌握地震波通过的地壳介质的波速信息,而这些信息很难确定,往往用小范围内的测量值代替整个路径上的地震波波速;二是计算震源到接收点的理论格林函数用时过长,虽然计算机内存不断在扩大,计算速度也不断在提高,但计算理论格林函数的用时依然不能满足需求。用经验格林函数方法模拟地震动的优势是使用具有与大地震相同的震源机制、相似传播路径的小震反演为格林函数,免去了计算理论格林函数的巨大计算量,也使得模拟的波形和反应谱都与地震记录符合得较好(药晓东等,2015李启成,2016),是1种值得推广的、有效的地震动模拟方法。随机方法模拟地震动所用的参数容易获得,计算速度快,但模拟结果针对的是大量地震动的平均效果,对具体的地震动记录的模拟不够理想,有些会有较大差异(王海云,2004)。

目前经验格林函数方法也有值得探讨的地方。首先是由于大地震的断层滑动分布一般是不均匀的,而实际模拟时很少考虑滑动位错不均匀的因素。虽然很多学者把凹凸体概念引进经验格林函数以减小模拟的误差(罗奇峰,1989Irikura等,1994Irikura,2006),但目前凹凸体概念与经验格林函数方法结合仅仅是用应力降描述凹凸体(Irikura等,1994, 2004),仍然没有考虑到大地震断层位错不均匀的因素。

其次,经验格林函数的小震数目必须是整数。但小震的数目是根据大小地震断层的长度、宽度和位错比值决定的,该比值一般不为整数,故实际应用时该整数的取值只能靠经验判断,而根据经验判断使得经验格林函数方法理论上不严密。

最后,经验格林函数方法强调大小地震必须满足“相似条件”(Irikura等,1994, 2004),即大小地震断层面长、宽和滑动位错是等比例的,这明显与实际矛盾,也限制了经验格林函数方法的适用范围。

针对目前经验格林函数方法存在的上述问题,用改进的经验格林函数方法(李启成,2010)与未改进的经验格林函数方法(Irikura,1983)对比模拟唐山大地震,以检验改进方法的有效性。

1 经验格林函数的改进

改进的经验格林函数方法计算公式如下(李启成,2010李启成等,2016):

$ {U_{{\rm{CB}}}}(x, t) = \sum\limits_{{\rm{l}} = 0}^{{N_L}} {\sum\limits_{m = 1}^{{N_W}} {\sum\limits_{k = 1}^{{N_D}} {\frac{{{R_{{\rm{CB}}}}(\theta, \phi)}}{{{R_{{\rm{CS}}}}(\theta, \phi)}} \times \frac{{{r_0}}}{{{r_{{\rm{lm}}}}}}} } } \times {K_{\rm{l}}} \times {K_{\rm{m}}} \times {K_{\rm{d}}} \times {V_{\rm{t}}} \times {U_{{\rm{CS}}}}(x, t - {t_{{\rm{sml}}}} - {t_{{\rm{dklm}}}}) $ (1)

其中,${U_{{\rm{CB}}}}(x, t)$是大地震加速度、速度或位移随时间的变化;${U_{{\rm{CS}}}}(x, t - {t_{{\rm{sml}}}} - {t_{{\rm{dklm}}}})$是小地震加速度、速度或位移记录KlKmKd分别为大小地震断层面长度比、宽度比和位错比的比值。

$ {N_L} = \left[ {\frac{L}{{{L_{\rm{S}}}}}} \right] $ (2)
$ {N_W} = \left[ {\frac{W}{{{W_{\rm{S}}}}}} \right] $ (3)
$ {N_D} = \left[ {\frac{D}{{{D_{\rm{S}}}}}} \right] $ (4)
$ {t_{{\rm{dklm}}}} = {r_{{\rm{lm}}}}/{v_{\rm{c}}} + \sqrt {{{(\xi - {\xi _0})}^2} + {{(\eta - {\eta _0})}^2}} /{v_{\rm{r}}} + ({N_D} - 1){\tau _{\rm{s}}} $ (5)
$ {t_{{\rm{slm}}}} = ({r_{{\rm{lm}}}} - {r_0})/{v_{\rm{c}}} $ (6)

LWD分别为大地震矩形断层面长、宽和位错;LSWSDS分别为小地震矩形断层面长、宽和位错;vc为剪切波速;vr为断层破裂速度,为剪切波速的0.8倍;${R_{{\rm{CB}}}}(\theta, \phi)$为大震的辐射因子;${R_{{\rm{SC}}}}(\theta, \phi)$为小震断层的辐射因子;${\tau _{\rm{s}}}$为小地震的上升时间;r0为小震震中到接收点的距离;rlm是大地震子断层到接收点的距离;tslm是由于大地震子断层和小地震断层到接收点的距离不同引起的时间延迟;tdklm是由于大地震断层破裂产生的时间延迟;NL是大小地震断层面的长度比;NW是大小地震断层面的宽度比;ND是大小地震的滑动位错比;Vt是大小地震断层的滑动速度比,该参数考虑了大小地震滑动速度不同的因素;ξ是大地震断层面沿走向方向的坐标,ξ0是沿走向的坐标原点;η是大地震断层面沿倾向方向的坐标,η0是沿倾向的坐标原点。

在公式(1)中,Vt参数的引入是考虑到大小地震断层的破裂速度不同的因素,使得模拟结果更符合实际情况。

一般情况下,大小地震长、宽和位错的比不为整数,此时KlKmKd取值为:当l=1,2……NL-1时,Kl=1;当1=NL时,Kl=大小地震断层长度比值的小数部分;当m=1,2……NW-1时,Km=1;当m=NW时,Km=大小地震断层宽度比值的小数部分;当d=1,2……ND-1时,Kd=1;当d=NDKd=大小地震断层位错比值的小数部分。

公式(1)中l取1到NL;m取1到NW;d取1到NDKlKmKd不仅可以取整数,也可以取小数,这样就突破了目前经验格林函数法必须满足大小地震相似的条件;KlKmKd可以取值不同,也就是考虑了大地震断层位错滑动不均匀特性;KlKmKd取值的规范化使得经验格林函数法更具理论意义。

2 用改进的经验格林函数方法模拟唐山地震地震动

20世纪70年代,中国唐山地区发生了MS 7.8大地震。地震震中位于39.63°N,118.18°E,震中烈度高达Ⅺ度,震源深度11km(国家地震局《一九七六年唐山地震》编辑组,1982)。通过选取与唐山大地震有相同震源机制、传播相似路径且信噪比高的小震记录作为经验格林函数,分别用改进的和未改进的经验格林函数方法进行模拟,模拟地震动时程及反应谱分别与实际地震动记录对比,讨论改进方法的有效性以及断层位错不均匀对地震动的影响。

2.1 唐山地震断层模型

唐山大地震主震断层主要由2个断层组成,其中1个位于南段,走向为北东30°,另1个位于北段,走向为北东50°,2个断层都主要具有右旋走滑特征(James,2002吴迪,2008)。唐山大地震主震断层的2个断层面总长度为114km,宽约20km,总地震矩为1.4×1027dyne·cm,平均应力降约为20bar,断层平均位错为138cm。唐山地区P波的传播速度为5.6km/s,S波为3.3km/s,通常,断层破裂速度取剪切波速的0.8倍,约为2.5km/s(周蕙兰,1985),整个断层的上升时间为0.6s(谢小碧,1988)。唐山主震断层的参数汇总于表 1

表 1 唐山大地震断层参数 Table 1 Fault parameters of the Tangshan earthquake

北京饭店和官厅水库的强震记录仪记录到了唐山地震主震(MS 7.8)和宁河强余震(MS 6.9)。选取MS 6.9余震在不同地点记录的加速度时程,作为经验格林函数。

根据宁河地震(MS 6.9)震源机制解和有关反演资料(张之立等,1980朱传镇,1985),选择宁河地震参数如表 2所示,其上升时间为0.22s。

表 2 宁河地震(MS 6.9)参数 Table 2 The Ninghe earthquake parameters (MS 6.9)

根据公式(1)计算得到主震各子断层所含子事件个数,如表 3所示。其中大小地震的长度比为3.16,取NL=4;经式(1)计算得到大小地震的宽度比值为2.2,取NW=3;ND的取值与NW相同。

表 3 大地震子断层含子事件个数 Table 3 The mumber of sub-events in the main earthquake

根据吴迪(2008)确定了2个凹凸体的位置和尺寸,如图 1所示。凹凸体1的尺寸为24km×8km,凹凸体2的9km×5km。凹凸体平均位错为2.76m,非凹凸体平均位错为1.22m。根据文献(李启成,2010),位于凹凸体1中2个子断层的Vt取1.28,凹凸体2中2个子断层的Vt取0.81,其余子断层均不是凹凸体,Vt为0.68。


图 1 唐山主震断层划分和凹凸体位置 Fig. 1 The fault division and the position of the convex and concave body of the Tangshan earthquak
2.2 模拟结果

放置在北京饭店一楼的强震仪(39.907°N,116.404°E)记录了唐山主震和宁河强余震,用宁河强余震记录作为经验格林函数,其东西向加速度记录如图 2所示。图 3为实际地震动记录、用未改进的经验格林函数方法模拟(模拟1)和用改进的经验格林函数方法模拟(模拟2)得到的唐山大地震东西向加速度基岩反应谱(反应谱的阻尼比为5%)。图 4为唐山大地震地震动记录时程、用未改进的经验格林函数方法模拟和用改进的经验格林函数模拟得到东西向地震动时程之间比较,可以看出:相对而言,用改进的经验格林函数方法模拟得到的基岩加速度反应谱与实际地震动记录反应谱符合得更好,但加速度时程改进效果不明显。


图 2 在北京饭店一楼记录的宁河地震东西向加速度 Fig. 2 The west-east acceleration of the Ninghe earthquake recorded at the first floor of the Beijing Hotel

图 3 3种方法获得的唐山大地震东西向反应谱比较 Fig. 3 Comparison of the response spectra in west-east direction of the Tangshan earthquake from three methods

图 4 3种途径获得的唐山大地震东西向加速度时程 Fig. 4 The west-east acceleration time histories of the Tangshan earthquake from the three approaches

放置在北京饭店一楼的强震仪记录了唐山大地震宁河强余震南北向地震加速度时程,如图 5所示,用该记录作为经验格林函数。图 6为在北京饭店一楼记录到的用未改进的经验格林函数方法和改进的经验格林函数方法模拟得到的唐山大地震南北向加速度基岩反应谱,图 7为用未改进的经验格林函数方法模拟得到和用改进的经验格林函数方法的唐山大地震南北向加速度时程。可以看出:相对未改进的模拟方法,用改进的经验格林函数方法模拟得到的基岩加速度反应谱与记录符合得更好。


图 5 北京饭店一楼记录到的宁河地震南北向加速度时程 Fig. 5 the north-south acceleration record of the Ninghe earthquake received at the first floor of the Beijing Hotel

图 6 3种方法获得的南北向反应谱比较 Fig. 6 Comparison of the response spectra in the north-south direction of the Tangshan earthquake from three methods

图 7 3种途径获得的唐山大地震南北向加速度时程 Fig. 7 Acceleration time histories of the Tangshan earthquake in the direction of north-south recorded from the three approaches

放置于北京饭店地下室的强震仪记录了唐山主震和宁河强余震记录,其中宁河强余震东西向地震记录如图 8所示,以该记录为经验格林函数。图 9为在北京饭店地下室记录到的用未改进的经验格林函数方法和用改进的经验格林函数方法模拟得到的唐山地震东西向加速度基岩反应谱。图 10为用未改进的经验格林函数方法和用改进的经验格林函数方法模拟得到的唐山地震东西向加速度时程比较。可以看出:相对而言,用改进的经验格林函数方法模拟得到的基岩加速度反应谱与记录符合得更好。


图 8 北京饭店地下室记录的宁河地震东西向加速度记录 Fig. 8 West-east acceleration time history of the Ninghe earthquake recorded at the basement of the Beijing Hotel

图 9 3种途径获得的唐山地震东西向反应谱比较 Fig. 9 Comparison of the response spectra of the Tangshan earthquake in direction of west-east from three approaches

图 10 3种途径获得的唐山地震东西向加速度时程 Fig. 10 The acceleration time histories of the Tangshan earthquake in the direction of west-east from the three approaches

放置在北京饭店地下室的强震仪记录到唐山地震宁河强余震南北向地震记录,如图 11所示,以该记录为经验格林函数。如图 12为在北京饭店地下室记录到的用经验格林函数方法和用改进的经验格林函数方法模拟得到的南北向基岩加速度反应谱,图 13为记录到的用未改进的经验格林函数方法和用改进的经验格林函数方法模拟得到的南北向加速度时程,可以看出,用改进的经验格林函数方法模拟得到的加速度反应谱相对用未改进模拟方法得到的结果,更接近地震记录。


图 11 北京饭店地下室记录的宁河地震南北向加速度记录 Fig. 11 The north-south acceleration of the Ninghe earthquake recorded in the basement of the Beijing Hotel

图 12 3种途径获得的唐山地震南北向反应谱比较 Fig. 12 Comparison of the response spectra of theTangshan earthquake in the direction of north-south at the Beijing Hotel basement from three approaches

图 13 3种途径获得的唐山地震南北向加速度时程 Fig. 13 The acceleration time histories of the Tangshan earthquake in the direction of north-south from the three approaches

放置在官厅水库的强震仪(39.467°N,118.208°E)记录了唐山主震和宁河强余震,记录到的宁河强余震南北向的加速度时程如图 14所示,以此为经验格林函数。图 15为用未改进经验格林函数方法和用改进的经验格林函数方法模拟得到的唐山地震南北向基岩加速度反应谱,图 16为相应的加速度时程。可以看出:用改进的经验格林函数方法模拟得到的基岩加速度反应谱与记录符合得更好。


图 14 官厅水库记录的宁河地震南北向加速度记录 Fig. 14 The north-south acceleration time history of the Ninghe earthquake recorded at the Guanting Reservoir

图 15 3种途径获得的唐山地震南北向反应谱比较 Fig. 15 Comparison of the response spectra of the Tangshan earthquake in direction of South-North from three approaches

图 16 3种途径获得的唐山地震南北方向加速度时程 Fig. 16 The north-south acceleration time histories of the Tangshan earthquake recorded from the three approaches

放置在官厅水库的强震仪记录到的宁河强余震竖直向加速度时程,如图 17所示,以该地震记录作为经验格林函数。图 18为用未改进经验格林函数方法和用改进经验格林函数方法模拟得到的唐山地震竖直向基岩加速度反应谱,图 19为获得的地震动时程。可以看出:相对于目前的模拟方法,用改进的经验格林函数方法模拟得到的基岩加速度反应谱与记录符合得比较好。


图 17 官厅水库记录的宁河地震竖直向加速度时程 Fig. 17 The acceleration time history of the Ninghe earthquake recorded at the Guanting Reservoir in vertical direction

图 18 3种途径获得的唐山地震竖直向反应谱比较 Fig. 18 Comparison of response spectra of the Tangshan earth quake vertical direction from the three approaches

图 19 3种途径获得的唐山地震竖方向加速度时程 Fig. 19 The vertical acceleration time histories of the Tangshan earthquake from the three approaches

从以上模拟结果可以看出,无论是目前的经验格林函数方法,还是改进的经验格林函数方法,2种方法反应谱模拟的结果相对于地震动时程模拟的结果都比较接近地震记录,用改进的方法模拟结果更接近实际地震动记录。为定量研究模拟结果与实际记录的差异,定义一定周期下模拟反应谱的幅值与记录反应谱的比值为模拟误差,简称误差。图 20为用未改进模拟方法计算得到的反应谱(实线)以及用改进的模拟方法计算的到的反应谱的平均模拟误差(虚线)。可以看出,在周期0—10s,改进的模拟方法计算得到的误差都小于用未改进的模拟方法计算得到的误差,用改进方法最大模拟误差不超过10%,而用未改进方法计算得到的误差最大值达到12%。


图 20 模拟误差 Fig. 20 Simulation error
3 讨论与结论

Irikura(1983)提出用经验格林函数模拟地震动,假定大地震的断层面的滑动位错是均匀的。近年来,随着大地震记录数据的增多,人们发现大地震断层滑动位错不均匀的现象是普遍存在的。如果依然用均匀的位错模拟地震动,无法反应真实的断层滑动。本文用真实的不均匀的断层面位错代替整个断层面均匀的滑动位错,使模拟结果更符合实际。

Irikura(1983)还假定大小地震具有相似性,也就是大小地震断层长、宽和位错的比值为同一数值,但这显然不具备普遍性。本文分别计算大小地震断层长、宽和位错比值,真实反映了大小地震断层长、宽和滑动位错的数值。

经验格林函数模拟地震动与随机方法模拟地震动相比,随机方法模拟的结果是众多的地震动的平均效果(王海云,2004)。对于具体的地震动模拟,随机方法有时会出现较大偏差(Beresnev等,1998),尤其体现在地震动持时上,由于随机方法地震动持时采用的是震源持时,模拟结果的持时一般小于记录持时,从这一点看经验格林函数方法模拟地震动具有明显优点。对于地震动的模拟效果,人们往往更注意反应谱或傅里叶幅值谱的模拟效果。

分别用未改进的经验格林函数方法和改进的经验格林函数方法模拟了唐山大地震多个观测点的地震动时程和反应谱,地震动时程模拟结果均与记录符合得不够理想。分析其原因,一是与模拟的理论有关,经验格林函数方法要求作为经验格林函数的小震应与大地震的震源机制相同,并且大小地震地震波通过的介质应相同。但实际上,上述2种条件绝对的满足是不现实的;二是唐山大地震断层复杂,其参数的误差必然是客观存在的,以上情况使得模拟地震动时程不够理想。

分别计算了用2种方法模拟地震动得到的加速度反应谱误差,在0—10s的周期内,即工程关注的周期内,用改进方法模拟的平均误差均比用未改进的方法模拟的小,用改进方法最大模拟误差不超过10%,而用未改进方法计算得到的误差最大值达到12%。说明改进的方法是更有效的模拟地震动方法。

致谢: 感谢审稿专家的修改建议。
参考文献
国家地震局《一九七六年唐山地震》编辑组, 1982. 一九七六年唐山地震[M]. 北京: 地震出版社.
李启成, 2010. 经验格林函数方法模拟地震动研究[M]. 哈尔滨: 中国地震局工程力学研究所.
李启成, 孙颍川, 闫晓丹, 2016. 用格林函数方法模拟卢龙地震[J]. 地球物理学进展, 31(5): 1947-1952.
罗奇峰, 1989. 近场加速度的半经验合成[M]. 哈尔滨: 国家地震局工程力学研究所.
王海云, 2004. 近场强地震动预测的有限断层震源模型[J]. 哈尔滨:中国地震局工程力学研究所.
吴迪, 2008. 基于凹凸体模型的地震动半经验合成研究[M]. 上海: 同济大学.
谢小碧, 1988. 理论地震图方法及其在研究震源过程方面的应用[M]. 北京: 中国科学院地球物理研究所.
药晓东, 章文波, 于湘伟, 2015. 2008年汶川8.0级大地震近场强地面运动的模拟[J]. 地球物理学报, 58(3): 886-903.
张之立, 李钦祖, 谷继成, 等, 1980. 唐山地震的破裂过程及其力学分析[J]. 地震学报, 2(2): 111-129.
周蕙兰, 1985. 浅源走滑大震震源过程的某些特征[J]. 地球物理学报, 28(6): 579-587. DOI:10.3321/j.issn:0001-5733.1985.06.004
朱传镇, 1985.唐山及其邻区的地震活动.见: 刘恢先主编, 唐山大地震震害.北京: 地震出版社, 1-16.
Beresnev I. A., Atkinson G. M., 1998. Stochastic finite-fault modeling of ground motions from the 1994 northridge, California, earthquake. Ⅰ. Validation on rock sites[J]. Bulletin of the Seismological Society of America, 88(6): 1392-1401.
Boore D. M., 2003. Simulation of ground motion using the stochastic method[J]. Pure and Applied Geophysics, 160(3-4): 635-676.
Irikura K., 1983. Semi-empirical estimation of strong ground motions during large earthquakes[J]. Bull Disas Prev Res Inst, Kyoto University, 33: 63-104.
Irikura K., 1986. Prediction of strong acceleration motion using empirical Green's function. In: Proceedings of the 7th Japan Earthquake Engineering Symposium. December 8, Tokyo, 151-156.
Irikura K., 1989. Prediction of strong ground motion using empirical Green's function-some problems of synthetic procedure of ground motion[J]. Disaster prevention research institute annuals, Kyoto University, 32(B-1): 41-52.
Irikura K., Kamae K., 1994. Estimation of strong ground motion in broad-frequency band based on a seismic source scaling model and an empirical Green's function technique[J]. Annals of Geophysics, 37(6): 1721-1743.
Irikura K., 2000. Prediction of strong motions from future earthquakes caused by active faults-case of the Osaka basin. In: Proceedings of the 12th World Conference on Earthquake Engineering. Auckland: IAEE.
Irikura K., Miyake H., Iwata T., 2004. Recipe for predicting strong ground motions from future large earthquakes[J]. Disaster prevention research institute annuals, Kyoto University, 47(A): 25-45.
Irikura K., 2006. Predicting strong ground motions with a "recipe"[J]. Bull disas prev res inst Univ Tokyo, 81: 341-352.
James M., 2002. Earthquake Source Studies Related to Strong Ground Motion Estimates[J]. Report of Japan's Ministry of Culture, Science and Technology: 664-670.