High-efficiency Fragility Analysis Method of NPP Equipment Based on Adaptive ANN
-
摘要: 传统的结构或设备易损性分析方法需要提供大量的数值模拟样本，这对于规模庞大的核电结构并不适用。为此，研究核电设备高效易损性分析方法，首先，基于拉丁立方法构建随机地震动-土-结构样本，采用高效时域SSI分区并行计算方法得到部分样本模型的地震响应；然后，采用有限的数值模拟结果训练人工神经网络模型（ANN），通过量化ANN预测误差和精确度指标，采用自适应算法进行后续数值模拟和ANN训练，直至满足精确度阈值要求。该方法可以优化计算样本的选择，控制数值模拟的样本数量，提高易损性分析的计算效率。此外，将ANN不确定性整合到易损性曲线计算公式中，分别基于对数正态假定的回归法和蒙特卡洛（MC）增量法对某核电设备进行了易损性分析，并验证了ANN不确定性量化方法的正确性。Abstract: Traditional methods of structural or equipment fragility analysis require a large number of numerical simulation samples, which is not applicable to large scale nuclear power structures. Therefore, an efficient technique for the fragility analysis of nuclear power equipment is developed in this paper. Firstly, the random ground motion-soil-structure samples are constructed based on the Latin method, and an efficient time-domain SSI partitioned parallel calculation method is used to obtain the seismic response of partial samples. Then, the neural network model (ANN) is trained with limited numerical simulation results, and by quantifying the prediction error and accuracy index of the ANN. An adaptive algorithm is used for subsequent numerical simulation and ANN training until the accuracy threshold requirements are all met. This technique can optimize the selection of calculation samples, control the number of samples for numerical simulation, and improve the calculation efficiency of fragility analysis. In addition, this paper integrates the ANN uncertainty into the calculation formula of the fragility curve, and conducts the fragility analysis of a nuclear power equipment based on the logarithmic regression method and the Mento Carlo incremental method, which verifies the ANN uncertainty.
-
Key words:
- Fragility analysis /
- Artificial neural network /
- MC method /
- SSI /
- Nuclear power equipment
-
引言
为实现“碳达峰，碳中和”的双碳目标，我国正大力发展核电、水电、风电等新型清洁能源，而我国地处两大活跃地震带，地震多发，如何保障地震环境下的能源工程安全问题，尤其是核电工程安全问题，已成为设计环节的重点。经过多年发展，概率地震风险评价（SPRA）已成为评价核电工程结构地震安全性的主要分析方法之一，结构和设备易损性分析是其中至关重要的一环。易损性分析方法最早可追溯到20世纪70年代，Kennedy等（1980，1984）基于现有确定性设计安全裕度提出了安全系数法，该方法从可获得的数据以及核电厂房设计信息出发，给出了较为成熟的核电厂结构及设备地震易损性评估程序，将核电结构及设备抗震能力表述为抗震能力中值和几个不确定性系数的乘积，但在工程应用中很大程度上依赖于经验数据和专家判断。随着计算能力的不断提高，基于数值模拟的易损性分析方法逐渐发展了起来，一般通过非线性时程分析获得样本数据，代表性方法有回归分析法、最大似然估计法以及增量动力分析法。其中，回归分析法（Cornell等，2002；Zentner等，2011）基于对数线性回归模型，适用于任何类型的数据，不需要对加速度时程进行缩放直至结构或设备失效，可以节省数值模拟成本，但同时若没有观察到结构或设备的失效意味着要进行外推，存在一定的异方差性。最大似然估计法（Shinozuka等，2000；Zentner，2010）尽管也不需要对加速度时程进行调幅缩放，但需要结构或设备失效之前和之后的样本数据，才能收敛并产生有意义的结果，该方法最适合于损伤模式为非好即坏的结构或设备，对于连续损伤模式的设备会导致信息缺失。增量动力分析法（IDA）（Vamvatsikos等，2002，2004）根据不确定性因素（如结构材料属性等）抽取部分样本并通过调幅一系列地震动直至结构破坏或设备失效确定易损性曲线，该方法需要重复大量数值模拟，相对来说较为耗时，但是最直接、最精确的分析方法，目前广泛应用于一般工程结构的易损性评估中（王中阳等，2018；杨贝贝等，2020）。
在评估核电厂地震危险性的过程中，美国规范ASCE/SEI 4—16 《Seismic analysis of safety-related nuclear structures》指出，应该从概率的角度进行评估，同时当核电厂地基不是硬基岩时，需要考虑土-结构相互作用。这需要同时考虑结构和土层材料特性的不确定性，大大增加了非线性时程分析的样本数量，使得本就耗时的易损性分析过程更复杂。近年来，涌现出一系列有效提高计算效率的易损性分析方法，通过建立元模型表达地震输入与结构输出之间的统计关系，代替传统物理模型以提高计算效率。Unnikrishnan等（2013）基于高维模型表示的响应面法，结合蒙特卡洛（Monte Carlo，以下简称“MC”）模拟生成易损性曲线，显著减少了数值分析次数。Calabrese等（2013）将人工神经网络（Artificial Neural Network，以下简称为“ANN”）应用于砌体码头结构地震易损性分析中，在修正好的ANN模型上直接执行MC过程，明显较在物理模型上快，并发现液化对重力式码头地震响应有明显的不利影响。Gehl等（2016）将贝叶斯网络应用于具有多风险荷载的桥梁系统4种失效模式的概率预测，验证了只要采用合适的建模策略，贝叶斯网络在求解大型复杂结构系统时具有足够的高效性和鲁棒性。Mangalathu等（2018）利用Lasso回归模型确定了对桥梁构件抗震需求影响最大的不确定变量。这些方法均旨在避开了重复的数值模拟过程，但大多数方法未对地震动特征代表值进行系统性的选择，同时缺乏对元模型预测值不确定性的量化评估。Wang等（2018a，2018b）基于ANN建立了元模型，对地震动特征选择以及元模型的不确定性量化进行了相关研究，但在易损性曲线中仅凸显了地震动的随机性，未考虑土层材料的不确定性。刘鸿泉等（2022）考虑土-结构相互作用效应，基于ANN对某核电厂设备进行了地震易损性分析，但未对数值模拟样本数量及选择给出指导意见。
本文构建了快速高效的核电设备易损性分析方法，包括采用拉丁立方法进行采样，控制总体样本数量；采用高效的时域SSI分区并行算法进行核电工程结构地震响应分析，提高单个样本数值模拟的计算效率；采用自适应的ANN训练方法，通过ANN模型预测误差的反馈，优化数值模拟的样本选择，控制了数值模拟的样本数量，最后基于优化的ANN模型执行MC模拟，利用增量法对核电设备进行高效的易损性分析。本文将上述方法应用于某核电结构电气设备，给出该方法的具体实现过程，并将该方法与传统对数回归法进行比较，验证了该方法的可行性。
1. 基于自适应ANN的设备易损性分析方法
易损性曲线估计了结构或设备在特定地震水平下的条件失效概率，对工程结构有重要的指导意义。在实际应用中，通常利用地震强度指标描述地震强度水平（如地表峰值加速度PGA），若将实际发生的地震动用某一重要的地震动特征值（IM）表示，则易损性曲线可表示为在此地震水平下，结构或设备地震工程需求参数（EDP）超越其失效阈值的条件概率为：
$$ {\text{P}}(\alpha ) = P(y > {y_{\rm{c}}}|\alpha ) $$ (1) $$ y = f(a(t),{{\boldsymbol{v}}_{\text{m}}}) $$ (2) 式中，
$\alpha $ 为某一地震动特征值表征下的地震动水平；$y$ 为地震需求值（如结构层间位移或某一位置的应力等）；${y_{\rm{c}}}$ 为结构或设备某一失效模式下的失效阈值；$a(t)$ 为地震动输入；${{\boldsymbol{v}}_{\rm{m}}}$ 为结构和土层材料属性中可变因素。在某一地震强度水平下，强震记录具有加速度峰值、持时、波动形式等方面的随机性，土层及结构的材料属性也具有一定程度的不确定性，其对地震需求值有不可忽视的影响，如何量化这些可变因素的不确定性成为了易损性分析的关键所在。增量分析方法采用MC模拟统计各种可变因素的不确定性，对结构输入与输出关系不作分布假设，是最简单直接的方法，其基本思想是在相同的地震动水平
$\alpha $ 下选取N条地震动进行结构地震响应分析，考虑地震动的随机性和材料的不确定性，统计出结构或设备的条件失效概率。由于核电工程结构较复杂、地震波时程较长等原因，导致精确的时域数值模拟过于耗时，直接采用MC法进行数值模拟代价大，因此本文利用ANN提出了高效的地震易损性分析方法，具体流程如图1所示。
图 1 基于自适应ANN的易损性分析方法流程Figure 1. Flowchart of fragility analysis method based on adaptive ANN初步选取部分样本进行数值模拟，利用数值模拟结果构建ANN模型代替物理有限元模型，同时采用ANN的自适应过程优化数值模拟的试验设计，为后续数值模拟的样本选择提供指导性意见，基于优化后的ANN模型结合MC法进行易损性分析。其中，初步样本选择采用拉丁立方采样法进行，拉丁立方采样本质上是分层抽样，其采样结果在一定程度上保证了样本结构与总体结构的相似性，因此可以在不增加样本数量的同时保证精度，是ASCE/SEI 43—05 《Seismic design criteria for structures, systems, and components in nuclear facilities》 推荐的方法。
1.1 SSI时域数值模拟流程
准确高效的数值模拟方法是易损性分析的基础，而出于计算效率的考量，考虑SSI效应的大型核电结构不适合采用传统的商业软件进行数值模拟，因此有必要开发高效的时域有限元分析方法。本文采用陈少林等（2010，2017）提出的时域分区并行算法进行SSI数值模拟，可大大提高数值模拟的计算效率，计算流程如图2所示，主要包括：
（1）选取基岩表面实测强震运动。选取的地震波应该与场地危险性分析得到的一致风险谱保持一致。
（2）采用等效线性化方法（ELM）考虑土层的非线性效应。本文采用SHAKE 91（Idriss等，1992）进行一维土柱的等效线性化，迭代生成与地震动强度水平相匹配的土层模型。
（3）基于等效线性化后的土层参数利用传递矩阵法计算场地自由场，作为后续SSI分析的输入。
（4）采用高效时域分区算法进行核电结构土-结构相互作用分析。
（5）对SSI分析结果进行处理得到目标地震需求参数。
1.2 ANN模型训练及不确定性量化
1.2.1 ANN模型训练
经过数值模拟得到样本数据后，即可利用样本数据训练ANN模型。BP神经网络主要结构如图3所示，主要包含输入层、隐藏层和输出层，训练网络的过程主要依靠反向传播算法最小化误差函数，以获得最优的网络系数（包括权重系数和偏置系数）。在网络训练前，需要确定ANN的超参数（如ANN层数、隐藏层的神经元数目等），其对网络逼近物理模型的能力具有至关重要的影响，通常需要设计者结合经验多次试验确定。通过众多学者检验，结构输入与输出关系在对数空间内并不具有强非线性，因此一般采用3层神经网络（Wang等，2018a）。而对于隐藏层神经元数目，笔者建议采用交叉检验的方法进行多次验证，选择泛化性能最优的隐藏层单元数。此外，作为ANN模型的输入变量，地震动特征值的选择对于易损性曲线的建立有较大影响，同时也影响ANN模型性能。选取的地震动特征值应尽量与地震需求参数相关，足以表达结构输入与输出之间的关系。地震动特征值可以通过输入与输出之间的半偏相关系数（姚俊，2010）以及交叉检验进行选择。在实际应用中，往往并不将全部样本数据均用于ANN训练过程，而是将可用的数据分为三部分以避免网络过度拟合，主要包括用于训练ANN的训练集，用于监视训练过程、避免出现过拟合的验证集和用于测试网络性能的测试集。
1.2.2 ANN模型不确定性的量化
有限的训练样本集并不一定能得到最优的ANN模型，其预测值和真实值之间存在一定偏差，量化ANN模型的预测误差对易损性分析至关重要，不仅是评价网络性能的重要依据，也是形成易损性曲线的重要因素。另外，在自适应算法中，预测误差也是网络优化以及样本选择的重要指标。
将ANN模型用函数
$f$ 表示，假设网络的最优权重系数和偏置系数为$ {{\boldsymbol{w}}^*} $ ，而经过训练得到的权重系数及偏置系数为$ {\boldsymbol{w}} $ ，当训练模型和最优模型足够接近时，可以通过泰勒展开逼近最优模型预测值：$$ f({x_i};{{\boldsymbol{w}}^*}) = f({x_i};{\boldsymbol{w}}) + {\boldsymbol{h}}_i^{\rm{T}}({\boldsymbol{w}} - {{\boldsymbol{w}}^*}) $$ (3) $$ {\boldsymbol{h}}_i^{\rm{T}} = \left(\frac{{\partial f({x_i},{\boldsymbol{w}})}}{{\partial {w_1}}},...,\frac{{\partial f({x_i},{\boldsymbol{w}})}}{{\partial {w_p}}}\right) $$ (4) 式中，
${{\boldsymbol{h}}_i}$ 为ANN预测值关于权重系数及偏置系数的梯度向量。那么，对于该结构系统的某一输入值$ {\boldsymbol{x}} $ ，其输出观测值$ {\boldsymbol{y}} $ 可利用ANN模型预测值进行表示：$$ {y_i} = f({x_i},{{\boldsymbol{w}}^*}) + {\varepsilon _i} = f({x_i},{\boldsymbol{w}}) + {\boldsymbol{h}}_i^{\rm{T}}({\boldsymbol{\hat w}} - {{\boldsymbol{w}}^*}) + {\varepsilon _i} $$ (5) 式中，
$\varepsilon $ 为最优模型与观测值之间的误差，${\boldsymbol{h}}_i^{\rm{T}}({{\hat {\boldsymbol w}}} - {{\boldsymbol{w}}^*})$ 为最优模型与训练模型之间的误差。经过论证（Chryssolouris等，1996；Zio，2006），上述2个误差相互独立且均服从标准正态分布：$$ {\varepsilon _i}{\text{～}}N(0,\sigma _{\rm{r}}^2) $$ (6) $$ {\boldsymbol{h}}_i^{\rm{T}}({{\hat {\boldsymbol w}}} - {{\boldsymbol{w}}^*}){\text{～}}N(0,\sigma _{\rm{u}}^2) $$ (7) $$ \sigma _{\rm{r}}^2 = \frac{1}{{N - p}}\sum\limits_{i = 1}^N {{{({y_i} - f({x_i},{\boldsymbol{\hat w}}))}^2}} $$ (8) $$ \sigma _{\rm{u}}^2 = \sigma _{\rm{r}}^2{\boldsymbol{h}}_i^{\rm{T}}{({{\boldsymbol{J}}^{\rm{T}}}{\boldsymbol{J}})^{ - 1}}{{\boldsymbol{h}}_i} $$ (9) $$ {\boldsymbol{J}} = [{{\boldsymbol{h}}_1},{{\boldsymbol{h}}_2},...,{{\boldsymbol{h}}_i},...,{{\boldsymbol{h}}_N}] $$ (10) 式中，
${\boldsymbol{J}}$ 为$p \times N$ 的矩阵，是整个训练集的雅可比矩阵；$p$ 为ANN权重系数及偏置系数总数；$N$ 为训练集样本总数。在SPRA方法中，用随机变量表征结构或设备的抗震能力，其中随机变量的不确定性又分为随机不确定性和认知不确定性，随机不确定性来源于表征随机事件的根本原因，如当指定地震强度水平时，地震记录的变异性以及结构材料的变异性，其不能通过增加数据而消除。认知不确定性来源于对随机变量的认识程度，可以通过增加模拟样本的数量减小，在易损性分析中可以生成一簇带有置信度的曲线。
而经过对ANN模型预测误差的量化分析，可以发现ANN预测误差的不确定性同样分为两部分，其中
${\sigma _r}$ 是由最优模型与真实模型误差引起的，不能通过增加训练样本数量而减小，主要来源于事件本身随机性质，如样本中地震动的随机性以及结构材料属性的随机性，属于随机不确定性；而$ {\sigma _u} $ 是由最优模型与训练模型之间的误差引起的，可以通过增加样本数量，不断优化模型而减小，主要来源于主观意识上对模型不确定变量认知的不足，属于认知不确定性。因此，可以通过ANN模型预测误差中的随机不确定性计算易损性曲线，而通过其中的认知不确定性确定曲线的置信区间，进一步可以通过自适应算法，减小ANN预测误差中的认知不确定性，提高易损性曲线的置信度。
1.3 ANN模型的自适应算法
训练后的ANN模型与最优模型之间存在误差，主要是因为样本数据不足，导致ANN欠拟合，较为直观的感受是在训练数据较稀疏处ANN预测结果会引起较大的误差，存在较大的认知不确定性。为提高ANN模型的泛化性能，即降低认知不确定性，需要增加用于训练的样本数据，采用自适应的训练算法，可以在少量的数值模拟样本基础上，判断哪些样本的认知不确定性最大，进而选择更合理的样本进行数值模拟，在一定程度上优化了数值模拟的试验设计，具体步骤如下：
（1）在全部地震动-土-结构样本集N中随机选择N0（N0<N）条样本作为初始样本集，进行数值模拟。
（2）利用数值模拟结果训练ANN模型。
（3）利用训练后的ANN模型计算全部样本集，对于任意ANN输入，可通过其预测误差计算每个预测值的误差指标
$ {\delta _i} $ ：$$ {\delta _i} = \left| {{s_i} - {s_{{\text{mean}}}}} \right|/{s_{{\text{mean}}}} $$ (11) $$ {s_i} = \sqrt {\sigma _{\rm{r}}^2 + \sigma _{{\rm{u}},i}^2} $$ (12) 式中，
${s_{{\text{mean}}}}$ 为整个训练集的平均预测误差。（4）根据误差阈值判断ANN模型预测值的精确度，定义误差阈值为训练样本集中最大的误差指标：
$$ {\delta _{{\text{crit}}}} = \max ({\delta _i}) $$ (13) 如果全部样本预测结果误差均小于误差阈值，则认为训练的ANN模型可以模拟全部样本，且精度在可控范围内，否则选择误差指标较大的若干样本进行数值模拟，并入训练样本集中，直至全部样本满足精度要求为止。
1.4 基于ANN计算易损性曲线
经过训练及优化后，可以用ANN模型代替物理模型进行易损性分析。相较于传统物理模型，ANN模型计算效率得到了提升，可以满足MC模拟的计算量需求，因此，不仅可以利用回归法计算易损性曲线，而且可以利用增量法计算易损性曲线。
采用回归法进行分析时，计算流程和传统的对数回归分析方法一致，不同的是数据是由ANN预测得到的而非数值模拟得到的，由于ANN预测值本身存在随机不确定性，因此计算曲线时需要同时考虑回归残差和ANN预测残差，其计算公式为：
$$ {P_f}(\alpha ) = \Phi \left(\frac{{c\ln \alpha + \ln b - \ln {y_{{c}}}}}{{\sqrt {\beta _{}^{\text{2}}{\text{ + }}\sigma _{\rm{r}}^{\text{2}}} }}\right) $$ (14) 式中，
$b,c$ 为回归系数；$\beta $ 为回归残差标准差；${\sigma _{\rm{r}}}$ 为ANN预测误差中代表随机不确定因素的标准差。采用增量法进行分析时，需要在相同的地震动水平
$\alpha $ 下执行MC模拟，选取N条地震动进行结构地震响应分析，考虑地震动的随机性和材料的不确定性，统计结构或设备的条件失效概率。其中，MC模拟可以利用ANN模型快速实现，仅需要利用地震动特征值的分布类型及边际范围，随机生成一系列特征值作为输入即可。同样的，在统计过程中，需要考虑ANN预测误差，其包括两部分，具有随机不确定性的误差，用于描述系统本身存在的随机性，如材料不确定性等；具有认知不确定性的误差，由于对模型缺乏认知导致的，用于生成一簇带有置信度的易损性曲线，计算公式如下：$$ {{\text{P}}_f}(\alpha ) = \frac{1}{N}\sum\limits_{i = 1}^N {b_i} $$ (15) $$ {{{\tilde P}}_f}(\alpha ) = \frac{1}{N}\sum\limits _{i = 1}^N {c_i} $$ (16) 式中，
${{y}_{c}} < {{\hat{y}}_{i}}(\alpha )+{{\varepsilon }_{{\rm{r}}}}$ 时，bi=1，反之为0；${{y}_{c}} < {{y}_{i}}(\alpha )+{{\varepsilon }_{{\rm{r}}}}+{{\varepsilon }_{{\rm{u}}}}$ 时，ci=1，反之为0。$\varepsilon_{\rm{r}} $ 为具有随机不确定性质的误差；${\varepsilon _{\rm{u}}}$ 为具有认知不确定性质的误差。当忽略认知不确定性时，式（16）可退化为式（15），代表中值易损性曲线。2. 实例分析
2.1 地震动-土-结构样本建立
本文参考某核岛结构模型（李小军等，2019）进行数值模拟。为了考虑土-结构相互作用，根据核电模型规模，建立了尺寸为640 m×360 m×60 m的场地模型，地基参数基于美国若干核电厂厂址调查假定的6种AP1000核岛设计地基中的上限软-中土场地（李小军等，2015），本文将该设计地基简化为3层水平场地，每层厚度均为20 m。为了满足模拟精度要求，场地模型离散为2 m×2 m×2 m的单元，在土层场地上部埋置刚性基础用于承载整个核电结构，刚性基础尺寸为92 m×60 m×16 m，总体模型如图4所示。其中，红色五角星为要进行易损性分析的目标电气设备，设备被简化为楼板上的附加质量，无法直接获取地震反应，因此通过楼板加速度反应谱量化设备反应，选取电气设备阻尼比为0.05，由于内部设备频率一般集中在5～33 Hz，因此本文选取目标设备所在节点处5～33 Hz的平均楼板加速度反应谱（AFSA）作为地震需求参数。
由于缺乏试验资料，电气设备失效阈值的获取是难点，本文着眼于高效设备易损性分析流程的建立，因此简单地以一种假定方法求解设备的失效阈值（尚昆，2014），该方法假定核电内部设备在设计地震动强度下有95%保证率不破坏，通过对地震数值模拟结果的概率统计确定设备的失效阈值。具体步骤如下：（1）依据核电结构的设计地震反应谱（DRS）选取多条强震运动，并根据结构基本周期处的设计谱值对地震动进行调幅；（2）对核电结构进行地震数值模拟（不考虑SSI效应），并获取相应的地震需求参数；（3）根据概率统计数值模拟结果，以全部数据的95分位值作为设备的极限能力值。本文通过对30条符合条件的地震动进行模拟和统计，并经过正态检验，发现结构的响应近似服从对数正态分布，经过统计分析，本文目标设备的极限能力值用AFSA衡量标定为0.5 g。
从概率角度进行设备易损性分析，需要考虑材料属性的不确定性以及地震动的随机性，本文在PEER网站上选取150条实测地震记录进行模拟，为使足够数量的地震动达到设备失效水平，对实测地震记录按照相同的方法进行调幅，调幅后的地震动反应谱如图5所示。考虑材料属性的不确定性，与随机地震动相对应，需要生成150个土-结构样本。采用有界的拉丁立方法进行抽样，考虑结构弹性模量以及土层剪切波速的变异性进行抽样建模。假设材料的分布类型为对数正态分布（Liel等，2009），根据其变异系数对各随机变量在95%置信区间内进行拉丁立方采样（即随机变量分布范围为
$u \pm 3\sigma $ ，其中$u$ 为样本均值，$\sigma $ 为样本方差），生成150个土-结构样本，并与随机地震动组合形成150个地震动-土-结构样本，结构及土层材料属性如表1所示。
图 5 设计反应谱及150条地震反应谱Figure 5. Design response spectrum and 150 seismic response spectrums表 1 核电结构及土层材料不确定性Table 1. Uncertainties in material parameters of NPP and soil结构 土层 结构类型 弹性模量/GPa 变异系数 名称 剪切波速/（m·s−2） 变异系数 NAB 24.7 0.2 L1 560 0.2 NSB 32.9 0.2 L2 673 0.2 SCV 210.0 0.2 L3 794 0.2 2.2 样本数值模拟
由于模型庞大，计算过程涉及近200万个单元的求解，本文采用陈少林等（2010，2017）提出的显-隐式结合方法进行求解，其中核电结构反应通过商业软件ANSYS求解，场地模型通过自编的FORTRAN程序进行求解，二者通过刚性基础进行耦合，具体原理及算法在文献中有详细介绍。
本文初步选取30条地震动-土-结构样本进行数值模拟。其中土层等效线性化采用李小军等（2015）提出的
${\text{G}}/{G_0} - \gamma $ 和$\lambda - \gamma $ 非线性特性曲线进行，并假设核电模型中的材料本构为线性本构，则在等效线性化后整个系统变成了线性系统，通过传递函数法减少不必要的计算耗时，以提高计算效率。2.3 神经网络训练及验证
针对30个地震动-土-结构样本数值模拟结果，采用MATLAB神经网络工具包建立ANN模型，模型结构为3层BP神经网络，其训练过程包括地震特征值的选择以及隐藏层单元数目的确定。选择地震特征值时需要区分主要特征值和次要特征值，主要特征值用于构建与地震需求参数的函数关系，表征地震强度水平，一般选取与地震需求值最相关的变量；其他次要特征值用于表达地震动的随机特性，在MC模拟中进行采样。更多的特征值选择更能还原地震动的随机特性，但要注意过多的特征值选择可能会降低网络的计算效率和泛化能力。本文推荐先根据相关性分析确定主要特征值，再利用交叉验证选择其他比较合适的次要特征值。
初步选取表2中的地震特征值作为候选，并通过相关性分析计算其与需求参数之间的相关系数及半偏相关系数，其中平均谱加速度（Average Spectral Acceleration，ASA）定义为阻尼比为0.05时地震动在5～33 Hz的平均加速度反应谱，其与DM的半偏相关系数最大，选择为主要的地震动特征值m0，其他影响较大的地震动特征值有PGA、Tp，理论上和DM有一定的相关关系，应在ANN模型建立过程中考虑其影响。根据计算结果构建以下4种工况作为输入进行交叉检验：（1）Case1：以ASA作为输入；（2）Case2：以ASA、PGA作为输入；（3）Case3：以ASA、PGA、Tp作为输入；（4）Case4：全部特征值作为输入。
表 2 地震动特征值Table 2. Seismic intensity measures地震动特征值 定义 R RSP 地震动特征值 定义 R RSP 峰值加速度（PGA） $ \max \left| {a(t)} \right| $ 0.73 0.18 累计绝对速度（CAV） $\displaystyle \int_0^{{t_{{\text{total}}}}} {\left| {a(t)} \right|{\rm{d}}t} $ 0.38 0.03 峰值速度（PGV） $ \max \left| {v(t)} \right| $ 0.29 −0.13 最大反应谱（PSa, max） $ \max (PS_{\rm{a}}(T)) $ 0.58 −0.03 峰值位移（PGD） $ \max \left| {d(t)} \right| $ 0.36 0.07 显著周期（Tp） $ \arg \max (PS_{\rm{a}}(T)) $ −0.26 0.16 阿里亚斯烈度（IA） $ \dfrac{{\text{π}} }{{2 g}}\displaystyle\int_0^{{t_{{\text{total}}}}} {a{{(t)}^2}{\rm{d}}t} $ 0.49 −0.03 平均谱加速度（ASA） $ \displaystyle\int_5^{33} {PS_{\rm{a}}(f){\rm{d}}f} $ 0.85 0.33 注：R为相关系数，RSP为半偏相关系数 交叉检验除了判断ANN的输入变量外，还可以用于确定隐藏层单元数目，本文采用10-折交叉检验进行模型验证，针对上述4种工况，分别以不同的隐藏层神经元数目进行训练，通过交叉检验判断模型的平均误差，选择最佳的ANN隐藏层单元数目。经过检验，无论是哪种工况，当h=1时网络模型最稳定，进一步地，当以Case3工况作为输入时，模型性能最稳定，如图6所示，因此本文选取隐藏层单元数目为1的3层BP神经网络进行训练，以Case3工况作为输入。
图 6 不同工况的10-折交叉检验（隐藏层单元数目为1）Figure 6. 10-fold cross validation for different cases （number of hidden layer neurons is 1）为验证对ANN模型预测误差量化的正确性，以Case3作为输入，将30个初始样本数据按照8∶1∶1的比例分为训练集、验证集、测试集，建立ANN模型。基于24个训练样本集，采用回归分析方法，分别基于FEM数据和ANN预测数据计算了设备的易损性曲线，如图7所示。由图7可知，基于ANN的曲线较基于FEM的曲线更陡峭，意味着其不确定性更小，主要原因是ANN训练过程中为防止过拟合，控制了数据的离散性。而当考虑了ANN预测误差中的随机不确定性时，按照式（14）将其修正到易损性曲线的计算公式中，得到的易损性曲线与FEM曲线完全重合，验证了ANN模型预测误差中随机不确定性因素量化方法的正确性。
2.4 ANN自适应过程及易损性曲线计算
ANN自适应过程避免了大量样本的数值模拟过程，本文在上述30个初选地震动-土-结构样本基础上进行ANN自适应过程，采用前文算法和步骤，选取每次增加的样本数目为5，经过反复迭代，直至全部样本的预测精确度满足条件为止，认为此时的ANN模型可以覆盖全部样本的输出预测。
利用训练后的ANN模型可以快速生成计算数据，适合利用增量法执行MC模拟，逐点计算失效概率，具体计算过程如下：
（1）主要特征值m0选为ASA，生成逐个计算点，根据设备抗震能力，在[1.8, 5.8]范围内线性抽取100个ASA特征值作为计算点，逐点计算失效概率。
（2）随机生成其他次要特征值IMs的组合，代表地震动的随机性。本文根据其他特征值的分布类型（对数正态分布）进行概率抽样，若其他特征值不止1个，则采用拉丁立方抽样，随机组合成1 000个特征样本，即针对每个IM0生成1 000个IMs样本计算设备的失效概率。由于未知特征值的边际分布，本文根据特征值的分布类型在其95%置信区间内进行采样，暂不考虑各特征值之间的相关性。
（3）逐点计算核电设备在IM0水平下的失效概率，根据式（8）、式（9）量化ANN模型的不确定性，结合式（15）、式（16）计算易损性曲线。
采用自适应算法，针对Case3工况建立ANN模型，其中迭代次数为2，即仅进行40次地震动-土-结构样本的数值模拟覆盖了全部样本的模拟任务。通过ANN自适应过程优化试验样本的选择，大大提高了计算效率。Case3工况下设备易损性曲线计算过程如图8所示。
由图8可知，随着迭代次数的增加，即随着样本数量的增加，ANN预测值最大误差指标越来越小，易损性曲线的置信区间逐渐变窄，意味着ANN模型预测数据的精确性得到提高，在样本数量增加不多的情况下可以较大程度地提高易损性曲线的可信度。同时发现，采用MC增量法计算的易损性曲线和采用回归法计算的易损性曲线有较好的相似性，而MC法为通过概率统计直接计算易损性曲线，因此该现象间接验证了回归法中对数线性回归模型假定的合理性。
总体而言，本文基于自适应ANN模型的核电设备高效易损性分析方法是可行的，基于数值模拟的设备易损性分析耗时主要来源于数值模拟过程，本文从控制数值模拟样本数量和提高数值模拟计算效率方面减少易损性分析耗时，其高效性体现在：（1）自适应ANN算法控制数值模拟样本数量，通过采用自适应ANN算法，将150个数值模拟样本降为40个数值模拟样本，整体分析效率提高了3.75倍；（2）采用分区并行算法进行SSI分析，其耗时相当于普通串行算法的26.6%（陈少林等，2017），通过提高单次样本数值模拟的计算效率，进一步将整体分析效率提高了近4倍。因此，综合分区并行SSI算法和自适应的ANN训练方法设备易损性分析效率是基于一般串行SSI算法的15倍左右。
3. 结论
本文建立了高效的设备易损性分析方法，并将其用于求解某核电设备的地震易损性曲线中。通过拉丁立方法建立了150个随机地震动-土-结构样本，采用高效时域SSI分析方法初步选取了30个样本进行地震响应分析，并通过自适应的算法训练ANN模型，最终基于40个数值模拟结果建立了ANN模型，并对核电设备进行了易损性分析，主要得出以下结论：
（1）ANN训练样本数据较稀疏处会带来较大的预测误差，主要来源于认知不确定性，该误差会使易损性曲线的置信区间变宽，而通过自适应的迭代过程增加样本数量，可以得到有效改善。
（2）对比了MC增量法和回归法计算的易损性曲线，从某种程度上验证了对数回归假定的合理性。
基于ANN模型的易损性分析方法充分利用了FEM模拟结果，尽可能地减少了计算样本数量，提高了计算效率，为发展重大工程结构设备易损性分析提供了参考。
-
图 1 基于自适应ANN的易损性分析方法流程
Figure 1. Flowchart of fragility analysis method based on adaptive ANN
图 5 设计反应谱及150条地震反应谱
Figure 5. Design response spectrum and 150 seismic response spectrums
图 6 不同工况的10-折交叉检验（隐藏层单元数目为1）
Figure 6. 10-fold cross validation for different cases （number of hidden layer neurons is 1）
表 1 核电结构及土层材料不确定性
Table 1. Uncertainties in material parameters of NPP and soil
结构 土层 结构类型 弹性模量/GPa 变异系数 名称 剪切波速/（m·s−2） 变异系数 NAB 24.7 0.2 L1 560 0.2 NSB 32.9 0.2 L2 673 0.2 SCV 210.0 0.2 L3 794 0.2 
下载: 导出CSV
表 2 地震动特征值
Table 2. Seismic intensity measures
地震动特征值 定义 R RSP 地震动特征值 定义 R RSP 峰值加速度（PGA） $ \max \left| {a(t)} \right| $ 0.73 0.18 累计绝对速度（CAV） $\displaystyle \int_0^{{t_{{\text{total}}}}} {\left| {a(t)} \right|{\rm{d}}t} $ 0.38 0.03 峰值速度（PGV） $ \max \left| {v(t)} \right| $ 0.29 −0.13 最大反应谱（PSa, max） $ \max (PS_{\rm{a}}(T)) $ 0.58 −0.03 峰值位移（PGD） $ \max \left| {d(t)} \right| $ 0.36 0.07 显著周期（Tp） $ \arg \max (PS_{\rm{a}}(T)) $ −0.26 0.16 阿里亚斯烈度（IA） $ \dfrac{{\text{π}} }{{2 g}}\displaystyle\int_0^{{t_{{\text{total}}}}} {a{{(t)}^2}{\rm{d}}t} $ 0.49 −0.03 平均谱加速度（ASA） $ \displaystyle\int_5^{33} {PS_{\rm{a}}(f){\rm{d}}f} $ 0.85 0.33 注：R为相关系数，RSP为半偏相关系数
下载: 导出CSV
-
陈少林, 唐敢, 刘启方等, 2010. 三维土-结构动力相互作用的一种时域直接分析方法. 地震工程与工程振动, 30（2）: 24—31Chen S. L. , Tang G. , Liu Q. F. , et al. , 2010. A direct time-domain method for analysis of three-dimensional soil-structure dynamic interaction. Journal of Earthquake Engineering and Engineering Vibration, 30(2): 24—31. (in Chinese) 陈少林, 王俊泉, 刘启方等, 2017. 基于显-隐式格式的三维时域土-结相互作用分析的异步并行算法. 中国科学: 技术科学, 47（12）: 1321—1330 doi: 10.1360/N092017-00212Chen S. L. , Wang J. Q. , Liu Q. F. , et al. , 2017. Asynchronous parallel algorithm for three-dimensional soil-structure interaction analysis based on explicit-implicit integration scheme. Scientia Sinica Technologica, 47(12): 1321—1330. (in Chinese) doi: 10.1360/N092017-00212 李小军, 侯春林, 戴志军等, 2015. 核岛结构设计地基场地及计算基底效应研究. 岩土力学, 36（8）: 2201—2208Li X. J. , Hou C. L. , Dai Z. J. , et al. , 2015. Research on site effects of soil layers and bedrock on designing the foundation of nuclear island structure. Rock and Soil Mechanics, 36(8): 2201—2208. (in Chinese) 李小军, 宋辰宁, 周国良等, 2019. 核岛结构PCS水箱FSI效应简化方法研究. 振动与冲击, 38（2）: 6—12, 32Li X. J. , Song C. N. , Zhou G. L. , et al. , 2019. Simplified method for simulating the FSI effect of PCS water tank in a nuclear island building. Journal of Vibration and Shock, 38(2): 6—12, 32. (in Chinese) 刘鸿泉, 陈少林, 孙晓颖等, 2022. 基于神经网络的核电厂设备易损性分析. 力学学报, 54（7）: 2059—2070Liu H. Q. , Chen S. L. , Sun X. Y. , et al. , 2022. Vulnerability analysis of NPP equipment based on neural network. Chinese Journal of Theoretical and Applied Mechanics, 54(7): 2059—2070. (in Chinese) 尚昆, 2014. 考虑SSI效应的核电厂安全壳及内部结构抗震能力评估. 哈尔滨: 哈尔滨工业大学, 15—17Shang K., 2014. Seismic assessment of nuclear power plants containment and internal structure considering SSI. Harbin: Harbin Institute of Technology, 15—17. (in Chinese) 王中阳, 车佳玲, 张尚荣等, 2018. 基于能量方法设计的RC框架结构易损性分析. 震灾防御技术, 13（3）: 524—533Wang Z. Y. , Che J. L. , Zhang S. R. , et al. , 2018. Seismic fragility analysis of RC frame structure based on energy balance. Technology for Earthquake Disaster Prevention, 13(3): 524—533. (in Chinese) 杨贝贝, 王志涛, 张秀彦, 2020. 基于IDA方法的框架结构震害风险评估. 震灾防御技术, 15（1）: 21—32Yang B. B. , Wang Z. T. , Zhang X. Y. , 2020. Seismic risk assessment of frame structures based on IDA method. Technology for Earthquake Disaster Prevention, 15(1): 21—32. (in Chinese) 姚俊, 2010. 半偏相关系数的计算公式及其应用. 常州工学院学报, 23（5）: 4—7 doi: 10.3969/j.issn.1671-0436.2010.05.002Yao J. , 2010. Calculating formula of semi-partial correlation coefficient and its application. Journal of Changzhou Institute of Technology, 23(5): 4—7. (in Chinese) doi: 10.3969/j.issn.1671-0436.2010.05.002 Calabrese A. , Lai C. G. , 2013. Fragility functions of blockwork wharves using artificial neural networks. Soil Dynamics and Earthquake Engineering, 52: 88—102. doi: 10.1016/j.soildyn.2013.05.002 Chryssolouris G. , Lee M. , Ramsey A. , 1996. Confidence interval prediction for neural network models. IEEE Transactions on Neural Networks, 7(1): 229—232. doi: 10.1109/72.478409 Cornell C. A. , Jalayer F. , Hamburger R. O. , et al. , 2002. Probabilistic basis for 2000 SAC federal emergency management agency steel moment frame guidelines. Journal of Structural Engineering, 128(4): 526—533. doi: 10.1061/(ASCE)0733-9445(2002)128:4(526) Gehl P. , D’Ayala D. , 2016. Development of Bayesian networks for the multi-hazard fragility assessment of bridge systems. Structural Safety, 60: 37—46. doi: 10.1016/j.strusafe.2016.01.006 Idriss I. M., Sun J. I., 1992. User's manual for SHAKE91: a computer program for conducting equivalent linear seismic response analyses of horizontally layered soil deposits. Davis: University of California. Kennedy R. P. , Cornell C. A. , Campbell R. D. , et al. , 1980. Probabilistic seismic safety study of an existing nuclear power plant. Nuclear Engineering and Design, 59(2): 315—338. doi: 10.1016/0029-5493(80)90203-4 Kennedy R. P. , Ravindra M. K. , 1984. Seismic fragilities for nuclear power plant risk studies. Nuclear Engineering and Design, 79(1): 47—68. doi: 10.1016/0029-5493(84)90188-2 Liel A. B. , Haselton C. B. , Deierlein G. G. , et al. , 2009. Incorporating modeling uncertainties in the assessment of seismic collapse risk of buildings. Structural Safety, 31(2): 197—211. doi: 10.1016/j.strusafe.2008.06.002 Mangalathu S. , Jeon J. S. , Desroches R. , 2018. Critical uncertainty parameters influencing seismic performance of bridges using lasso regression. Earthquake Engineering & Structural Dynamics, 47(3): 784—801. Shinozuka M. , Feng M. Q. , Lee J. , et al. , 2000. Statistical analysis of fragility curves. Journal of Engineering Mechanics, 126(12): 1224—1231. doi: 10.1061/(ASCE)0733-9399(2000)126:12(1224) Unnikrishnan V. U. , Prasad A. M. , Rao B. N. , 2013. Development of fragility curves using high-dimensional model representation. Earthquake Engineering & Structural Dynamics, 42(3): 419—430. Vamvatsikos D. , Cornell C. A. , 2002. Incremental dynamic analysis. Earthquake Engineering & Structural Dynamics, 31(3): 491—514. Vamvatsikos D. , Cornell C. A. , 2004. Applied incremental dynamic analysis. Earthquake Spectra, 20(2): 523—553. doi: 10.1193/1.1737737 Wang Z. Y. , Pedroni N. , Zentner I. , et al. , 2018 a. Seismic fragility analysis with artificial neural networks: application to nuclear power plant equipment. Engineering Structures, 162: 213—225. doi: 10.1016/j.engstruct.2018.02.024 Wang Z. Y. , Zentner I. , Zio E. , 2018 b. A Bayesian framework for estimating fragility curves based on seismic damage data and numerical simulations by adaptive neural networks. Nuclear Engineering and Design, 338: 232—246. doi: 10.1016/j.nucengdes.2018.08.016 Zentner I. , 2010. Numerical computation of fragility curves for NPP equipment. Nuclear Engineering and Design, 240(6): 1614—1621. doi: 10.1016/j.nucengdes.2010.02.030 Zentner I. , Humbert N. , Ravet S. , et al. , 2011. Numerical methods for seismic fragility analysis of structures and components in nuclear industry - application to a reactor coolant system. Georisk: Assessment and Management of Risk for Engineered Systems and Geohazards, 5(2): 99—109. doi: 10.1080/17499511003630512 Zio E. , 2006. A study of the bootstrap method for estimating the accuracy of artificial neural networks in predicting nuclear transient processes. IEEE Transactions on Nuclear Science, 53(3): 1460—1478. doi: 10.1109/TNS.2006.871662 -
下载:
点击查看大图
计量
- 文章访问数: 145
- HTML全文浏览量: 79
- PDF下载量: 14
- 被引次数: 0
