基于响应面和重要抽样法的隧道衬砌结构时变可靠度

姚贝贝，孙 钧，2

（1.同济大学 地下建筑与工程系，上海200092；2.杭州丰强土建工程研究院，浙江 杭州310006）

地下结构的安全可靠与否，不仅影响其正常运营使用，而且直接关系到人民的生命财产安全.而盾构隧道建造费用高，服役期长（一般为100年），破坏危险性大，因此，要绝对保证盾构隧道的安全可靠.在目前的研究中，抗力不随时间变化的可靠度研究多，抗力随时间变化的可靠度研究少.结构的耐久性不足时会造成结构抗力的降低，从而使结构的安全性下降[1].显然，在结构耐久性成为国际土木界所关注问题的今天，研究结构时变可靠度是有现实意义的.

结构时变可靠度是指将结构在设计使用期限内的抗力和作用效应看成是随机变量，结构抗力随时间不断降低的可靠度分析方法，它是一个相当复杂的问题.隧道衬砌结构受力大部分处于偏心状态下，极限状态函数的非线性程度较高，特别是在考虑抗力随时间变化的情况下.目前，国内外对结构可靠度的常用分析方法主要有一次二阶距法、Monte Carlo方法、随机有限元法、响应面法[2－3].一次二阶距方法运算简单，但主要是针对功能函数能够明确表达的结构，对于复杂结构而言常难以写出功能函数的显式形式；Monte Carlo方法其模拟的收敛速度与基本随机变量的维数无关，极限状态函数的复杂程度与模拟的过程无关，具有直接解决问题的能力.近几年，还发展了各种效率更高的抽样方法如对偶抽样，重要性抽样[4]、分层抽样、控制变数法等，但对于实际工程的结构破坏概率通常小于10－3以下量级的范畴时，该方法的模拟数目就会相当大，占据大量的计算时间，效率较低；随机有限元方法是另一种手段，但是需要对确定性结构分析程序加以改造，要形成一个通用的随机有限元程序来描述工程实际中各种随机性，目前尚有一定困难[5]；响应面法具有思路清晰，方法简便，计算量小的特点，是目前最有发展前景的结构可靠度分析方法之一，但该方法计算的是可靠度指标是几何可靠度而不能完全代表结构的真实可靠度[6].

由于本文研究的是结构时变可靠度，极限状态方程较为复杂，故将结合响应面法和重要抽样两种方法进行时变可靠度计算.即利用响应面法计算出设计验算点和响应面方程，再引入重要抽样法以得到的设计验算点作为抽样中心，以响应面方程作为近似功能函数进行抽样，并考虑结构抗力随时间衰减，从而求解隧道衬砌结构服务寿命期内的时变可靠度.

1 响应面法

响应面法起源于实验设计，而后用于结构可靠度的数值模拟.其基本思想就是对于隐含的或需花费大量时间确定的真实的功能函数或极限状态面，用一个容易处理的函数（称为响应面函数）或曲面（称为响应面）代替.一般为先设计一系列变量值X，每组变量值构成一个试验点即样本点xi，并逐步计算结构相应的一系列功能函数值Zi（i＝1，2，…，n），构造变量组和功能函数值之间的明确函数关系，利用其近似代替真实功能函数，再利用常用的方法计算结构可靠度.响应面法的关键是响应面函数对取样点的很好拟合，响应面函数通常选取二次多项式，本文选取的近似响应面函数为不含交叉项的非完全二次多项式：

式中：a0，bi，ci（i＝1，2，…，n）为待定系数.

确定响应面函数的关键是确定待定系数，先根据室内试验得到的各试样得到的变量均值作为初始迭代点，利用Matlab编程进行循环确定最接近均值点的一系列变量点，从而求出待定系数，得到一个时间变量对应一个响应面函数，以此反复进行迭代，直至两次迭代结果在误差范围之内.通常在设计试验点时采用二水平因子设计或中心复合设计，二水平因子设计是取因子的上水平和下水平，当有n个因子时，需要2n次试验.中心复合设计是在二水平设计的基础上，在增加原点和2n个坐标轴上的点.式（1）沿坐标轴代表真实功能函数，试验点可沿坐标轴在均值点μX附近选择，其中沿坐标轴Xi轴的试验点具有的坐标为xi＝μXi±fσXi，其中f＞0，是一任意因子.这是一种只有在坐标轴上的点的中心复合设计.在这2n＋1个点处计算原真实功能函数的值，并由此得到线性方程组，解之可得到2n＋1个未知系数a0，bi，ci.具体的求解步骤如下：

（1）假定初始迭代点x＝（x1，x2，…，xn）T，一般取平均值μX.

（2）选取f值，一般取1，2，3，本文程序中选取f＝2.

（3）通过结构数值分析或试验在各个展开点处计算功能函数的估计值（i＝1，2，…，2n＋1），并形成相应的系数矩阵，利用结构数值试验方法，通过Matlab编制专用程序进行计算.

（4）利用式（1）求解待定系数a0，bi，ci（i＝1，2，…，n）.

（5）计算时变可靠度指标β及验算点x＊.

（6）计算在验算点处的功能函数的估计值.

（7）通过线性插值公式可得到新的x，即利用公式

（8）返回步骤（3）进行迭代，直至前后两次β值相差＜ε，‖x＊‖＜ε.

2 重要抽样法

设结构的功能函数为Z＝gX（X），基本随机变量X的联合概率密度函数为fX（x），则结构的失效概率表示为

可靠度指标与失效概率的对应关系

式（2）—（3）中：I（x）为x的指示函数（或称特征函数、示性函数），规定当x＜0时I（x）＝1，反之，I（x）＝0；Df是与gX（X）相对应的失效区域；Φ（·）为标

准正态分布的累积概率函数.

利用Monte Carlo方法表示式（2）可写为

式中，N为抽样模拟总数.

抽样模拟总数N可近似的表示为

式（5）意味着抽样数目N和成反比，而工程结构中的失效概率通常是较小的，这说明N必须要有足够大的数目才能给出正确的估计，很明显，只有利用方差缩减技术，降低抽样模拟数目N，才能使Monte Carlo法在实际工程可靠性分析中得以应用.因此赵国藩等学者提出了效率更高的重要抽样方法.

直接的Monte Carlo抽样法得到的随机变量X的样本点xi（i＝1，2，…，N）多集中在联合概率密度函数fX（x）的最大值附近，该点一般比较接近X的均值点μX.实际的结构失效应为小概率事件，从而μX处于可靠域而不在极限状态面上，在失效域内的样本点很少，实现一次Z＜0的机会很小.因此，重要抽样方法的基本思想是通过改变随机抽样中心，使样本点有较多机会落入失效域，增加使功能函数Z＜0的机会.

假定存在一个重要抽样概率密度函数pV（v），满足下列关系

则式（2）可变为

重要抽样法就是选用pV（v）进行抽样，可能改变原抽样的重要区域，增加样本点落入失效域的机会，但若绝大部分落入失效域内也对求解不利.因此，可以将抽样中心取在失效域内对结构失效概率贡献最大的点v＊，即最可能失效点.v＊可通过以下最优化问题求解：

若用pV（v）对V抽样，得到样本vi＝（vi1，vi2，…，vin）T，i＝1，2，…，N，则pf的无偏估计值为

pV（v）的基本变量V的各量为正态随机变量，V的方差可取对应的原随机变量X的方差的1～2倍，V的均值可取最大可能点v＊或验算点x＊.

综上所述，重要抽样法就是先以 Monte Carlo直接抽样得到的样本点为初始抽样中心，利用式（8）进行最优化设计，得到最可能失效点v＊，将v＊重新作为抽样中心，进行迭代计算，利用式（7）可进行失效概率的计算.该过程均只是Matlab程序中实现.

根据可靠度指标与失效概率一一对应的关系，可有式（9）求出失效概率，从而求出可靠度指标.

3 基于响应面法与重要抽样法结合的可靠度计算方法

对于实际工程中，一般极限状态方程为隐式，直接利用重要抽样法，其抽样中心和抽样区域很难确定，并且直接采用重要抽样法，抽样次数仍然比较高，效率较低；而响应面法对于隐式方程能较精确地确定其验算点和验算点附近拟合精度较高的响应面.因此，将两者结合可提高效率和计算精度，以响应面法得到验算点作为重要抽样法的抽样中心，从而构造抽样函数进行重要抽样，计算结构失效概率和可靠度指标.具体步骤如下：

（1）假定初始迭代点x＝（x1，x2，…，xn）T，一般取平均值μX.

（2）选取f值，一般取1，2，3，本文程序中选取f＝2.

（3）通过结构数值分析或试验在各个展开点处计算功能函数的估计值i（i＝1，2，…，2n＋1），并形成相应的系数矩阵.

（4）利用式（1）求解待定系数a0，bi，ci（i＝1，2，…，n）.

（5）计算可靠度指标β及验算点x＊.

（6）若满足收敛条件｜β（k）－β（k－1）｜＜ε，则输出β，否则以为样本中心，返回步骤（3）进行迭代，直至满足收敛条件.

（7）输出验算点x＊，以x＊作为重要抽样中心和抽样函数pV（v）的均值，以原随机变量方差的1～2倍作为抽样函数的均值；

（8）对nV个随机变量选取nS个正态分布的随机抽样点Vi；

（9）计算随机抽样点对应的功能函gX（vi）；

（10）以原随机变量的均值和方差构造随机抽样点对应的联合概率密度函数fX（vi）；

（11）以验算点为均值，以原随机变量方差的一倍为重要抽样函数的方差构造重要抽样概率密度函数pV（vi），并计算fX（vi）／pV（vi）的数值；

（12）计算功能函数值gX（vi）＜0出现的总次数；

（13）计算失效概率pf.

通过以上原理和步骤，本文利用Matlab软件，结合隧道的工程实际概况，编制了服务寿命期限内，结构时效可靠度和可靠度指标的计算程序，并通过实例分析验证了其正确性.

4 实 例

某越江公路隧道全长8 955.26 m，隧道内设计车速80 km·h－1，结构设计使用年限为100年.其中江中段东线长7 471.65 m，西线长7 469.363 m，为双管盾构隧道.江中圆隧道的衬砌外径为15.0 m，内径13.7 m，环宽2.0 m，环厚0.65 m.工程区年平均温度15.7℃，年平均相对湿度0.8.

本文研究时效可靠度，因此要考虑抗力和作用响应随时间变化的统计特征，计算时取单位宽度进行计算，通过有限元计算，本文以衬砌结构偏心受压状态为例，利用近似概率法建立衬砌结构的极限状态方程.在结构服务寿命期限内荷载作用效应S（t）不能超该期限内结构的抗力R（t），即：

当不满足式（10）时结构将失效，因此结构的极限状态方程为

当Z（t）＞0时，结构安全；当Z（t）＜0时，结构失效.

4.1 抗力随时间衰减的随机计算模型

对隧道衬砌结构，其受力状态一般为偏心受力.根据文献[7]规定，公路隧道衬砌结构矩形截面的强度计算公式，并结合文献[8]推导的构钢筋混凝土矩形截面小偏心受压构件，部分截面受压时的抗力表达式，以期使计算结果更为准确.当构件小偏心受压时，文献[9]的推导的抗力表达式变为：

式中：Ra为混凝土的极限强度；h0为有效高度；Ag，Ag1分别为受拉区、受压区钢筋截面面积；e，e1分别为轴向力作用点到受拉钢筋、受压钢筋合力点的距离；Rg1为钢筋的抗压计算强度；这些参数假定为正态分布，统计特征按照文献[10]取值.另外三个参数为定值，取值分别为b＝1 000 mm，B＝－1 309 MPa，C＝1 047 MPa.

本文在计算钢筋锈蚀后衬砌结构强度时，考虑环境因素和构件老化导致的抗力随时间衰减，并引入损伤系数.分别反映钢筋锈蚀引起的钢筋截面损失及强度降低、混凝土截面损伤、钢筋和混凝土间粘结力下降三方面的损伤效应.钢筋锈蚀后截面强度计算则需考虑损伤系数as和ac的影响，受拉和受压区分别变为asRaAg和asacRaAg1.

在公式（11）的基础上，引入表1中的损伤系数，并考虑时间效应，得到隧道衬砌结构在服务寿命期限内，结构抗力随时间衰减的计算公式为：

式中：Ra（t）为混凝土的极限强度；as（t）为钢筋锈蚀后钢筋强度和截面面积降低系数随机过程；ac（t）为钢筋锈蚀后钢筋与混凝土的协同工作系数；这3个参数均与时间有关，通过室内氯离子渗透和扩散试验利用Fick第二定律得到的氯离子扩散系数，和室内钢筋电化学快速锈蚀试验得到的钢筋锈蚀量来确定，其相关参数参考文献[11]中的数值，均假定为正态分布，由于数值较多，所以在文中不在列出.h0c含义见表1.

表1 随机参数统计特征Tab.1 Statistical characteristics of random parameters

4.2 作用效应S统计特征

本文在进行荷载效应求解过程中将围岩的力学性能、衬砌材料的力学性能、衬砌结构的几何尺寸、计算模式等看作随机变量，其统计特征见表1.利用“荷载－结构”法，采用Ansys软件中的可靠度计算专用模块对衬砌结构进行内力分析，则可得到截面弯矩M0和轴力N0的均值和方差.

4.3 考虑抗力随时间衰减的极限状态方程

考虑抗力随时间衰减的极限状态方程如下：

通过以上极限状态方程，利用本文编制的程序，求解出该越江隧道衬砌结构从建成开始使用到服务寿命100年内的时效可靠度指标和对应的失效概率，运算结果见表2.

表2 服务寿命期内可靠度指标与失效概率Tab.2 Time－dependent reliability indexes and failure probability in the service lifetime

整理出衬砌结构服务寿命期内时变可靠度指标与服务时间的关系图（图1）.

图1 衬砌结构时效可靠度指标与时间关系图Fig.1 Time－dependent reliability indexes and time diagram of tunnel lining

从图1可以看出，在结构100年服务寿命内，结构时效可靠度随时间是不断衰减的，且与时间成指数关系，通过回归分析，表示如下：

在工程结构中，对于缺少实测资料的工程，在基于近似概率法进行耐久性设计时，求解耐久性设计系数η则是关键，而耐久性系数则是时效可靠度指标的函数；两者一一对应的关系如下：

式中：β0为结构初始可靠度指标，一般按文献 [12]取值，本文取3.7；βt结构达到服务寿命时的可靠度指标，一般取βt＝β0.

从式（15）可以看出，耐久性系数与时效可靠度指标是相互对应的关系，通过本文的方法程序，可计算出结构时效可靠度指标，进而求得耐久性系数，才可对拟建结构在缺少实测数据的情况下进行耐久性设计.本文通过式（16）计算出该隧道在达到服务寿命100年时，结构的耐久性系数为η＝0.733 9，结合近似概率法，可验证该隧道能够满足100年的服务寿命.可见，本文的方法程序，在结构进行耐久性设计时是有现实意义的.

5 结论

本文通过考虑隧道环境因素影响及时间效应，基于隧道衬砌结构抗力随时间的衰减，结合近似概率法得到的极限状态方程，利用响应面法和重要抽样法计算衬砌结构的时效可靠度，主要有以下结论：

（1）介绍在高度非线性极限状态方程情况下，将响应面和重要抽样法结合的理论计算方法，该方法为隐式非线性功能函数的越江隧道工程可靠度分析提供了一种有效的方法；

（2）在隧道服务寿命期内，考虑抗力随时间衰减，利用衬砌结构极限状态方程，通过室内快速试验得到随机参数，以近似概率法得到的极限状态方程为基础，通过响应面法得到的验算点和精度较高的响应面函数，以该结果作为重要抽样的抽样中心，利用Matlab软件编制了隧道在服务寿命期内时效可靠度求解的专用软件；并回归分析了隧道衬砌结构时效可靠度与时间呈指数衰减规律；

（3）通过时效可靠度指标与耐久性系数的一一对应关系，利用可靠度指标求解耐久性系数，为越江隧道衬砌结构的耐久性设计提供了基础.

[1] 赵国藩，金伟良，贡金鑫.结构可靠度理论[M].北京：中国建筑工业出版社，2000.ZHAO Guofan，JIN Weiliang，GONG Jinxin.Theory of structural reliability[M].Beijing：China Architecture and Building Press，2000.

[2] Faravelli L.A response surface approach for reliability analysis[J].Journal of Engineering Mechanics，1989，115（12）：2 763.

[3] Bucher C G，Bourgund U.A fast and efficient response surface approach for structural reliability problems[J].Structural Safety，1990，7（1）：57.

[4] Bucher C G，Bourgund U.A fast and efficient response surface approach for structural reliability problems[J].Structural Safety，1990，7（1）：57.

[5] 王建军，于长波，李其汉.工程中的随机有限元方法[J].应用力学学报，2009，26（2）：297.WANG Jianjun，YU Changbo，LI Qihan.Stochastic finite element methods in engineering[J].Chinese Journal of Applied Mechanics，2009，26（2）：297.

[6] 苏永华，方祖烈，高谦.用响应面方法分析特殊地下岩体空间的可靠性[J].岩石力学与工程学报，2000，19（1）：55.SU Yonghua，FANG Zulie，GAO Qian.Reliability analysis of special uderground rockmass by response surface method[J].Chinese Journal of Rock Mechanics and Engineering，2000，19（1）：55.

[7] JTG D70—2004公路隧道设计规范[S].北京：人民交通出版社，2004.JTG D70—2004 Code for design of road tunnel[S].Beijing ：China Communications Press，2004.

[8] 李敏.海底隧道衬砌结构可靠性研究[D].北京：北京交通大学建筑工程系，2005.LI MIN.The reliability of tunnel lining structure[D].Beijing：Department of Building Engineering of Beijing Jiaotong University，2005.

[9] Thomas M D A，Bamforth P B.Modelling chloride diffusion in concrete；ash and slag[J].Cement and Concrete Research，1999，29：487.

[10] GB50068—2001建筑结构可靠度设计统一标准[S].北京：中国建筑工业出版社，2001.GB50068—2001 Unified standard for reliability design of building structures[S].Beijing：China Architecture and Building Industry Press，2001.

[11] 陈海明.越江隧道衬砌结构耐久性设计若干关键技术研究[D].上海：同济大学地下建筑与工程系，2009.CHEN Haiming.Study on some key technologies on the durability design of tunnel lining structure under river[D].Shanghai：Department of Geotechnical Engineering of Tongji University，2009.

[12] GB50010—2002混凝土结构设计规范[S].北京：中国建筑工业出版社，2002.GB50010－2002 Code for design of concrete structures[S].Beijing：China Architecture and Building Industry Press，2002.