崔 臻，盛 谦，宋艳华，魏 倩

（1.中国科学院武汉岩土力学研究所 岩土力学与工程国家重点实验室，武汉 430071；2.电子科技大学 物理电子学院，成都 610054）

1 引 言

对于目前我国西南地区众多坐落在高山峡谷地带的高水头、大容量电站，地下厂房是一种经济的甚至是唯一的选择。相比于几何形态较为简单的隧道等地下工程，由于功能性的要求，大型水工地下洞室群通常空间关系较为复杂，同时由于西南地区活断层多、规模大、活动强度高，这些大型水利水电工程项目大部分处于西南强地震多发地区和高抗震设防地震烈度地区，设计地震动参数水平较高。这就使得这些大型地下洞室群的地震作用下的稳定性，成为工程设计中无法回避的关键科学问题。

目前针对大型地下洞室群地震响应的分析方法，多为采用若干条确定地震波进行时程分析[1－5]。在取得了大量的研究成果的同时，也应注意地震波是一种强烈的随机过程[6]，在数值计算中输入的具体地震波及其产生的地震响应均仅能视作随机过程的若干确定性表现，较难揭示其客观规律。

运用随机振动理论，将输入地震波和地震响应均视作随机过程，直接分析其统计规律，可以更好地反应地下洞室群地震响应的随机性特征和本质规律[7－9]。此外，随机地震响应分析也可以较方便与目前流行的概率设计理论相结合。因此，该方法对地下洞室群抗震设计理论的进一步完善和发展具有较显著的意义。

传统上岩土工程的随机地震动响应多采用脉冲响应函数进行分析[8－10]，本文试图直接将传递函数引入地震随机响应研究中。本文首先论述了采用传递函数计算地下洞室群地震响应随机参数的合理性，继而针对金沙江白鹤滩水电站地下厂房13#机组剖面，根据计算到的传递函数进行地震随机响应研究，并将随机响应计算结果与基于传递函数的确定性估算结果及实际时程结果进行了对比。最后，首次采用超越破坏理论，计算了白鹤滩水电站地下厂房洞室群在不同失效概率下对应的地震响应极值。上述研究可为地下洞室群的地震动力响应研究提供一些新的思路。

2 传递函数

传递函数源自于自动控制理论，是一种简洁有效的信号分析数学工具。传递函数反映了系统的输入与输出之间的频域关系，体现系统的固有特性，是系统在频域中的一个重要特征量。传递函数概念的引入，对振动理论起到了很大的推动作用[11－13]。

2.1 传递函数原理

对受动荷载 f （ t ）激励的简单质量(m)-弹簧(k)-阻尼(c)系统建立强迫振动动力学模型，其微分方程为

对式（1）两边进行拉普拉斯（Laplace）变换，并假设初值为0，得

式中：s 为拉普拉斯变换因子； ( )x s 为x 的拉普拉斯变换，相应 ( )f s 为f 的拉普拉斯变换。

其倒数称为机械导纳，又称为传递函数

传递函数在傅氏域中的变换 ( )H ω 称为频率响应函数，所以频响函数实际上是传递函数在傅氏域中的别称，在多数实际应用中对这2 个名词也不加以严格区别[13]。

2.2 地震响应频谱特性的传递函数分析

传递函数以频域形式表达了被研究系统的动力特性，从频域形式描述了系统（地下洞室群工程区域）对输入信号（地震波）的传递特性。地壳深部传来的地震波在岩体内部各处阻抗不同的区域存在着反射、折射和衍射等复杂的波动效应，使得地震波信号中的不同频率分量在岩体介质中的传播规律存在一定区别，一些频率分量被放大的同时，一些分量被抑制，使得输出信号相对于原输入信号的频谱特性有了新的分布特征。故传递函数在频域表述上的优势，使其可以用来研究地震动力响应的频谱特性。

2.3 确定性地震波输入下的响应估算

对于线性定常系统，若取得1 个已知传递函数，则系统的频域响应规律可以通过输入信号的频谱特性加以确定[14－16]。采用传递函数估算岩体对输入地震动的频域响应，将拟采用的地震波的傅立叶谱实/虚部与已知的传递函数实/虚部进行逐频复数相乘，即得到估算的傅立叶谱值，如下式：

式中： F (ω )为输入加速度的傅立叶谱；zj为数值模型中某点； F (ω , zj)为该点加速度时程的傅立叶谱估算值。

由于 估算 得 到的 傅立 叶 谱 F (ω , zj)中 的幅 值和相位信息是完整的，利用傅立叶逆变换方法IFFT[14－15]可得到时域响应 F ( t , zj)的估算值为

3 基于传递函数的地下洞室群随机地震响应

地震动通常分为上升段（初震），平稳段（强震）和下降段（衰减）3 个阶段，是明显的非平稳随机过程。对于地震动，工程中有2 种处理方法：一种为确定性方法，即采用某次实测地震加速度记录作为输入，计算系统的响应，但不能保证另一次地震能得到同样的结果；另一种为随机振动方法，即探讨地震随机过程的一般规律，一般而言，地震动在强震阶段的水平分量大致是平稳的。为简化分析，常常进一步简化为0 均值的平稳高斯随机过程。

随机振动的最主要特征是激励和响应事先不能用时间的确定函数表达，但仍可以利用统计的方法研究其规律性。相关研究的基本问题都是由输入的自相关函数或功率谱密度函数来确定体系输出的自相关函数或功率谱密度函数，从而确定体系响应的方差和均方差。

下文将证明，若已知1 个受随机激励作用下的线性定常系统的传递函数，其随机响应的均值、均方值、最大值等统计参数均可以用传递函数进行表达[6,11－12]。

3.1 随机响应的均值

仍以一简单质量-弹簧-阻尼系统为例，将激励更换为随机力 ( )F t ，动力学方程为

将 F (t )视作一系列脉冲激励的叠加，在t = τ至τ+ dτ 的 微 小 间 隔 内 激 励 产 生 的 脉 冲 冲 量 为F （τ ）dτ 。则式（7）的响应特性可用传递函数 H (ω )或其傅立叶逆变换脉冲响应函数 h （ t ）描述，而其杜哈姆积分形式解为

当 ( )F t 为平稳随机过程时，其稳态响应也是平稳随机过程。 ( )x t 的均值可表示为

将式（9）中的脉冲响应函数变换为传递函数，则响应的均值可表示为

可见响应的均值和激励的均值只相差一个常数。若激励均值为0，响应的均值也为0。

3.2 随机响应的均方值

根据式（8），并用λ1和λ2表示积分变量，随机响应的互相关系数可写为

由维纳-辛钦关系，平稳随机过程的自功率谱密度是其自相关函数的傅立叶变换，则随机响应的自功率谱密度表示为

注意方括号中即为激励的自谱 SF(ω )，且再次将脉冲响应函数转换为传递函数，有

式（13）表明，根据激励力的自谱和系统的传递函数即可求出响应的自谱。由自功率谱密度和随机过程均方值的关系（范马克关系）可求得随机过程的均方值：

3.3 随机响应的极值

均值为0 的平稳高斯随机过程 x ( t )在单位时间内，极值x = a或 r = a /σx的概率分布函数由达文波特关系决定为

式中： Ω = σx˙σx， σx˙为 x (t )导数过程的均方值，Td为激励持续时间。其概率密度函数为

由此可得随机响应极值的均值和均方值近似为

3.4 地下洞室群地震动力可靠度探讨

随机地震分析的最终目的是为了进行系统的动力可靠度分析。因此，在上文随机地震响应研究的基础上，对地下洞室群进行可靠度分析。首次超越破坏理论是目前结构动力可靠度分析中广泛采用的一种方法[8,10,17]，它是以随机响应过程 ( )x t 在规定时间内不超越限值x a= 的概率来定义系统的动力可靠度。

将地下洞室群的稳定性系数表示为广义结构抗力R，将地震荷载表示为广义荷载效应S，在二者相互独立且服从正态分布的假定下。功能函数Z= R S- 也服从正态分布。对此有可靠度指标β 和结构失效概率fP 分别满足

式中： Y = ( Z - μz)σz。常用的失效概率与可靠度的对应关系见表1。

表1 常用的失效概率与可靠度的对应关系 Table 1 Corresponding relationship between Failure probability and reliability

综上所述，较之传统上常用的利用脉冲函数进行随机地震分析的方法，基于传递函数进行地下洞室群系统随机地震响应时，仅需进行1 次计算，除可得到洞室群的随机地震响应以外，还可得到洞室群的地震响应频谱特性，并可以直接利用传递函数得到确定性地震动输入下洞室群的响应值用以进行验证，显示其对问题研究有更全面的优势。

4 工程实例

金沙江白鹤滩水电工程地下厂房为大跨度、高边墙的大型地下洞室群，主厂房跨度为32 m，工程位于地震活动强烈的高山峡谷地区，地震基本烈度Ⅶ度。选取右岸厂房13#机组剖面，进行随机地震响应研究。

数值模型、岩体参数等条件与文献[18]相同，研究中考虑高山峡谷地形、层间错动带、岩性和洞室群开挖等因素的影响。在3 大洞室各部位设置了18 个监测点，监测点分布和编号如图1 所示。

在以下的分析中，限于篇幅，仅给出了加速度指标的研究成果，但同理可以较容易的扩展到速度、位移及应力等指标。

图1 关键点布置及其编号 Fig.1 Key points layout and number

4.1 传递函数的计算及频谱分析

在不对输入地震波做更多假设的前提下，有限带宽随机白噪声是对地震动随机过程的一种较简单描述。白鹤滩水电站场址50 a 超越概率5%的加速度峰值为219 gal，由于本文中将岩体考虑为线性系统，故将峰值归一化为1 m/s2，以便于与峰值变化进行对比。按归一化剪切波输入有限带宽随机白噪声，得到洞室各部位的响应之后，算得传递函数。图2 为主厂房各监测点及尾调室层间错动带出露部位上下盘监测点的加速度绝对传递函数幅值。

图2 典型部位的加速度传递函数 Fig.2 Transfer function of acceleration for typical position

从图2 可见，地下洞室群各部位的传函幅值与常见的多自由度结构或水平场地的传递函数幅值相比，峰值较不明显，与地下结构不具备明显的自振频率的一般认识相符合，显示出地下结构更复杂的动力特性。

由传递函数幅值可进行洞室的地震响应频谱分析。结果表明，地震波自下而上向地表传播的过程中，不同频段的分量传播规律不尽相同，高频分量衰减较大。由于地震波的衍射作用，主厂房底板的响应总体上强于顶拱部位，底板部位以1.0～3.0 Hz为主要影响频段，以1.5～2.5 Hz最为强烈，在1.8 Hz附近具有1 个显著频率；顶拱部位以1.5～4.0 Hz为主要影响频段，以2.0～3.5 Hz最为强烈，在2.0 Hz及3.4 Hz 附近具有2 个显著频率，且在3.0～4.0 Hz频段，顶拱部位的响应强于底板部位，显示了洞室各部位不同的动力特性。这些结论基本与文献[18]采用小波包频谱分析手段得到的结论一致。

4.2 随机地震动响应

在地震随机响应和可靠度分析中，地震荷载常以功率谱密度表达。Kanai-Tajimi 模型是一种应用较广的地震动加速度平稳自功率谱，国内亦称金井清模型[6]。该模型表达式如下：

式中：gω 和gξ 分别为岩土体的特征频率和特征阻尼比。

作为比较，采用了基于Kanai-Tajimi 模型改进的平稳过滤有色噪声模型，即欧进萍模型[8]作为比较。欧进萍模型的表达式如下：

式中：hω 为反映岩土体特性的谱参数。

计算结果见表2、3。由表可见，在地下洞室群随机地震响应中，地震波自下而上传播，在洞室底板自由面发生反射现象，当下行反射波与上行入射波发生叠加时，地震波强度提升，表现为加速度峰值增加。而在地震波传播至洞室空腔之后，在洞室群上方“阴影区”具有较明显的衍射绕行现象，使得洞室顶拱附近区域地震波的幅值显示出一定程度的减少。在层间错动带附近，由于构成层间错动带的岩屑力学参数弱、阻抗低，地震波峰值同样有一定程度减小。这些结论与前文中基于传递函数频谱特性分析中得到的结论较为一致。

比较表2、3 中的结果可见，Kanai-Tajimi 模型的结果比欧进萍模型大7%左右。为了对结果进行验证，采用了本文2.3 节中提出的确定性地震估算方法，对洞室群在集集地震波和Kobe 地震波两条确定性地震波作用下的响应值进行了估算。同时采用实际时程计算对估算值本身进行了验证。如图3所示。

表2 Kanai-Tajimi 模型随机地震响应值 Table 2 Stochastic seismic response of Kanai-Tajimi model

表3 欧进萍模型随机地震响应值 Table 3 Stochastic seismic response of OU Jinping model

以主厂房底板为例，Kanai-Tajimi 模型得到的地震响应的可信极值为1.282±0.240 m/s2，欧进萍模型的为1.176±0.285 m/s2。以传递函数估算方法得到的确定性响应值分别为1.300、1.140 m/s2，而实际时程计算的结果分别为1.280、1.220 m/s2。以上结果表明，2 种地震动随机模型得到的结果均较为良好，Kanai-Tajimi 模型的结果相对保守一些。

值得提到的是，对2 条确定性地震波进行实际时程计算，在Dell-T5500 工作站上分别花费了40、20 d；若需要研究其他地震波，还需要重新进行大规模计算。而基于传递函数的估算和随机响应分析需要的仅是1 次时程计算和耗时几十秒的数学计算。理论上来讲，得到传递函数的时程计算可以非常短，但为了保证传递函数具有足够的FFT 变换点以维持精度，建议时程的采样点数最好多于1 024个，如此FFT 变换点可不少于2 048 个[19]。即当采用常用的0.02 s 采样频率时，计算传递函数的时程时间以大于21 s 为宜。

图3 估算及实际时程计算得到的加速度时程响应曲线 Fig.3 Acceleration response of estimated & calculated results

4.3 地震动力可靠度的探讨

以地下洞室群遭遇Ⅷ度烈度的地震为例，讨论地下洞室群的动力可靠度。假定50 a 设计基准期内场址最大地震烈度的概率分布服从极值Ⅲ型分布，场地地震烈度与地面地震加速度峰值之间存在如下关系[17]：

式中：mla 为地震烈度为I 的加速度峰值。由此可以得到洞室群各部位随机地震响应值与烈度之间的关系。图4 为主厂房底板随机地震响应值与场地地震烈度的关系。

图4 主厂房底板随机地震响应值与地震烈度的关系 Fig.4 Relationships between stochastic seismic response & seismic intensity of the main cavern floor

表4 场地烈度Ⅷ度时不同可靠度对应的极值 Table 4 Extremum of various reliability corresponding to seismic intensity of Ms 8

可见，随着可靠度指标增高，主厂房各部位响应的极值也随之增加。由表4 可见，例如需要保证主厂房顶拱在Ⅷ度烈度时失效概率小于1‰，则需要考虑的设计地震动加速度为2.9 m/s2；若需失效概率小于1‰，则需要考虑的设计地震动加速度为3.1 m/s2。这表明设计时若需要使得洞室群具有较高的抗震可靠性，则设计加速度指标需很高。过高的可靠度要求将会带来难以承受的工程成本，需要按照对危险的接受程度合理选择可靠度指标。

5 结 论

（1）较之传统常用的利用脉冲函数进行随机地震分析的方法，本文提出的基于传递函数研究地下洞室群系统随机地震响应的方法除可得到洞室群的随机地震响应以外，还可同时得到洞室群的地震响应频谱特性，并可以直接利用传递函数得到确定性地震动输入下洞室群的响应值进行验证，显示本文提出的方法对问题研究具有更加深刻的优势。

（2）本文方法得到的随机地震响应结果与基于传递函数的确定性估算结果与实际时程计算得到的确定性结果吻合良好。这说明较之计算成本高昂的确定性时程分析，本文采用的随机响应分析和基于传递函数的确定性估算结果的方法可以更加简明、有效地把握地下洞室群在地震作用下的动力响应整体规律。

（3）探讨了采用首次超越破坏理论对地下洞室群进行地震动力可靠度分析的可行性，为地下工程抗震可靠度计算中如何以动力方式考虑地震荷载的随机性提供了一些新思路。采用该理论，得到了白鹤滩水电站地下厂房洞室群在不同失效概率下对应的地震响应极值。

（4）传递函数、随机地震响应的理论目前仅建立在线性系统假定上，当需要考虑岩体的非线性力学特征时，其应用途径当作进一步研究。同时，如何将地震动的非平稳特性考虑进来，也是一个值得进一步研究的问题。

[1] 黄胜, 陈卫忠, 杨建平, 等. 地下工程地震动力响应及抗震研究[J]. 岩石力学与工程学报, 2009, 28(3): 483－490. HUANG Sheng, CHEN Wei-zhong, YANG Jian-ping, et al. Research on earthquake-induced dynamic responses and seismic measures for underground engineering[J]. Chinese Journal of Rock Mechanics and Engineering, 2009, 28(3): 483－490.

[2] 王如宾, 徐卫亚, 石崇, 等. 高地震烈度区岩体地下洞室动力响应分析[J]. 岩石力学与工程学报, 2009, 28(3): 568－575. WANG Ru-bin, XU Wei-ya, SHI Chong, et al. Dynamic response analysis of rock underground caverns in highly seismic region[J]. Chinese Journal of Rock Mechanics and Engineering, 2009, 28(3): 568－575.

[3] 李小军, 卢滔. 水电站地下厂房洞室群地震反应显式有限元分析[J]. 水力发电学报, 2009, 28(5): 41－46. LI Xiao-jun, LU Tao. Explicit finite element analysis of earthquake response for underground caverns of hydropower stations[J]. Journal of Hydroelectric Engi- neering, 2009, 28(5): 41－46.

[4] 赵宝友, 马震岳, 丁秀丽. 不同地震动输入方向下的大型地下岩体洞室群地震反应分析[J]. 岩石力学与工程学报, 2010, 29(增刊1): 3395－3402. ZHAO Bao-you, MA Zhen-yue, DING Xiu-li. Seismic response of a large underground rock cavern groups considering different incident angles of earthquake waves[J]. Chinese Journal of Rock Mechanics and Engineering, 2010, 29(Supp.1): 3395－3402.

[5] 张志国, 肖明, 陈俊涛. 大型地下洞室地震灾变过程三维动力有限元模拟[J]. 岩石力学与工程学报, 2011, 30(3): 509－523. ZHANG Zhi-guo, XIAO Ming, CHEN Jun-tao, et al. Simulation of earthquake disaster process of large-scale underground caverns using three-dimensional dynamic finite element method[J]. Chinese Journal of Rock Mechanics and Engineering, 2011, 30(3): 509－523.

[6] 胡聿贤. 地震工程学[M]. 北京: 地震出版社, 2006.

[7] 陈健云, 胡志强, 林皋. 超大型地下洞室群的随机地震响应分析[J]. 水利学报, 2002, (1): 71－75. CHEN Jian-yun, HU Zhi-qiang, LIN Gao. Random seismic response analysis of large-scale underground caverns[J]. Journal of Hydraulic Engineering, 2002, (1): 71－75.

[8] 严松宏. 地下结构随机地震响应分析及其动力可靠度研究[博士学位论文D]. 成都: 西南交通大学, 2003.

[9] 张鸿儒, 曹国安. 地下结构随机地震响应分析[C]//第五届全国土动力学学术会议论文集. 大连: 大连理工大学出版社, 1998.

[10] 李志强, 张鸿儒, 陆传波. 挡土结构物随机地震响应分析及动力可靠度研究[J]. 工程地质学报, 2006, 14(3): 345－350. LI Zhi-qiang, ZHANG Hong-ru, LU Chuan-bo. Stochastic analysis of seismic responses and dynamic reliability of retaining structures[J]. Journal of Engi- neering Geology, 2006, 14(3): 345－350.

[11] 季文美. 机械振动[M]. 北京: 科学出版社, 1985.

[12] 刘延柱. 振动力学[M]. 北京: 高等教育出版社, 1998.

[13] 傅志方. 模态分析理论与应用[M]. 上海: 上海交通大学出版社, 2000.

[14] BERND GIROD, RUDOLF RABENSTEIN, ALEXANDER STENGER. Signals and systems, 2nd ed[M]. New York: Wiley, 2001: 50.

[15] YAMAGUCHI Y, IWASHITA T. Earthquake motions at an embankment dam base and an estimation method of incident seismic waves using the observations[C]//The 32nd Joint Meeting of the U.S.-Japan Cooperative Program in Natural Resources Panel on Wind and Seismic Effects. Gaithersburg: [s. n.], 2001: 253－262.

[16] 蒋良潍, 姚令侃, 吴伟, 等. 传递函数分析在边坡振动台模型试验的应用探讨[J]. 岩土力学, 2010, 31(5): 1368－1374. JIANG Liang-wei, YAO Ling-kan, WU Wei, et al. Transfer function analysis of earthquake simulation shaking table model test of side slopes[J]. Rock and Soil Mechanics, 2010, 31(5): 1368－1374.

[17] 李桂清. 结构动力可靠性理论及其应用[M]. 北京: 地震出版社, 1993

[18] 崔臻, 盛谦, 刘加进, 等. 基于小波包的地下洞室群地震响应及其频谱特性研究[J]. 岩土力学, 2010, 31(12): 3525－3530. CUI Zhen, SHENG Qian, LIU Jia-jin, et al. Spectrum characteristics analysis of seismic response for under- ground houses via wavelet packet[J]. Rock and Soil Mechanics, 2010, 31(12): 3525－3530.

[19] The MathWorks Inc.. MATLAB signal processing toolbox user's guide[M]. Natick: The MathWorks Inc., 2010.