基于响应面法和麻雀搜索算法的结构有限元模型修正

徐 喆，辛景舟,2*，唐启智，肖维娜，李双江

(1.重庆交通大学省部共建山区桥梁及隧道工程国家重点实验室，重庆 400074；2.广西交通投资集团有限公司，南宁 530022；3.贵州毕节高速发展有限公司，毕节 551700)

建立精准的有限元模型是进行大型桥梁结构系统识别、损伤检测和承载力评估的基本条件[1]。然而，在设计阶段，采用确定性参数建立有限元仿真模型对实际工程结构设计建模时，需要对结构材料属性、几何特性以及边界条件进行一定的简化和近似处理，按照设计图纸建立的理想化模型很难完全准确地反映实际结构；此外，在施工阶段，工程建设人为误差、材料性能、荷载形式和环境条件的随机性等不确定性因素对结构性能具有显著影响[2], 这些不确定性也无法通过模型定量描述，导致了有限元模拟结果与结构表现的实际行为之间的误差超过工程允许范围。因此有必要根据实测结果对有限元模型进行修正，为结构后期运营期间的健康监测、管养维护、状态预测与评估提供科学理论指导。

目前，结构有限元模型修正方法按测试对象主要可划分为矩阵型修正法和参数型修正法两类[3]。矩阵型修正法由于在工程应用中的复杂性和低效性使得其应用具有一定的难度[4]。相反，参数型修正法通过灵敏度分析法[5-7]和响应面法[8-11](response surface model, RSM)调整设计参数来实现。对于实际工程结构，修正参数较多，采用灵敏度分析法需要大规模的反复迭代计算有限元模型，计算效率较低；而响应面法以简单低阶数学模型代替响应与各参数的物理系统关系，进而实现隐式关系的显式表达，在计算效率上更具有优越性。近年来，随着有限元模型修正研究成果不断涌现，基于响应面替代模型的结构有限元模型修正逐渐走进人们的视野，响应面法仍然是结构有限元模型修正中应用最为广泛的一种替代模型技术[12]，具有极高的发掘潜力和应用价值。

结构有限元模型修正是通过调整模型材料弹性模量、容重、截面刚度及边界条件等设计参数来实现的，使修正后的模型计算值与实际结构响应误差最小[13]。实际上，模型修正问题可转化为非线性优化问题，在设计参数约束空间内寻找一个全局最优解，其目标函数是由结构计算值与实际值之间的残差函数构成，优化后的目标函数越小，说明模型与实际结构状况越契合，修正效果越好。所以，为得到客观反映结构真实响应的有限元模型，选择合理的修正算法就显得尤为重要。传统的数学优化算法如梯度下降法、牛顿法、拉格朗日乘数法等，当目标函数维度较高、非线性较强时，算法容易陷入局部最优，寻优效果得不到保证。面对实际工程结构，响应与参数间往往表现出较强的非线性，应用上述方法进行有限元模型修正时，其精度无法得到保障。而智能仿生优化算法是模拟自然界生物群体的一般特性而形成自适应全局最优化随机搜索算法，可以有效避免传统算法的弊端，显著提高算法效率，近年来被广泛应用于结构有限元模型修正领域[14-17]。

麻雀搜索算法(sparrow search algorithm, SSA)是Xue等[18]于2020年提出的一种模拟麻雀群体觅食行为和反捕食行为的新型群体智能优化方法。该算法比较新颖，结构简单，具有全局寻优能力强，收敛速度快，迭代过程稳定，对高维数据不敏感等特点，一经推出，就表现出极高的应用前景[19]。碍于目前现有的常规优化方法在精度和效率已无法满足结构模型修正问题，有必要展开基于麻雀搜索算法在结构有限元模型修正方面的研究。

为寻求解决结构模型修正的新理论方法，现提出基于RSM-SSA的结构有限元模型修正方法，首次将麻雀搜索算法应用到有限元模型修正中。首先，依据结构构造特点及力学特性确定待修正参数和结构响应，构造结构响应与参数的二阶多项式响应面模型，依此联合动力响应计算值与实际值之间的残差函数；其次，通过由种群初始化、发现者和跟随者位置更新、侦查预警四部分组成的麻雀搜索算法，对目标函数进行迭代优化，在设计空间内寻找参数最优解；最后，通过悬臂梁数值算例，验证所提方法的可靠性与适用性。

1 基于RSM的目标函数构建

目标函数构建以参数尽可能地满足实际结构宏观响应需求为目标，将结构有限元模型修正机理由数学问题转化为优化问题，即对式(1)求解最优解：

(1)

式(1)中：d为待修正参数；Yc、Ye分别为实际结构响应和基于模型特征量响应；VLB、VUB为结构参数上下限；Y(d)为残差。

如何实现结构实际响应与参数复杂非解析隐性关系的显式表达，成为优化问题求解的关键。通过响应面法实现优化目标函数的构建，结合统计理论和数学建模技术建立输入与输出映射函数关系，以简单低阶的数学模型近似逼近结构实际响应特征量。

1.1 响应面函数的选择与试验设计

选用不含交叉项的二阶多项式作为响应面模型的基本函数形式，可表示为

(2)

样本点数量和分布特征显著影响着响应面模型构建效率及精度。样本数量太少，结构响应与参数之间的复杂隐性关系得不到完全映射；而样本数量较多，在一定程度上提高了拟合精度，但又从客观上延长了计算分析的时间，降低了试验效率。拉丁超立方设计(Latin hypercube design, LHD)很好地处理了拟合精度与试验成本这一矛盾问题，它能以较小的样本反映总体的变异规律，往往能有效改善样本的均值和方差，提高抽样效率和精度，抽样的次数可大大减少。为此，本试验设计选择基于多维分层抽样思想的拉丁超立方设计。

1.2 响应面的建立

将式(2)以矩阵形式表示为

Xλ=Y

(3)

式(3)中：

λ=(ab1b2…bnc1c2…cn)T；N为样本个数；Y=(Y1Y2…YN)T。

根据最小二乘原理，求得基函数的系数矩阵为

λ=(XTX)-1XTY

(4)

随着设计参数数量和多项式展开阶次的增加，响应面模型中的待定基函数系数急剧增加。然而，并非所有多项式展开项对响应的贡献都是显著的，因此如何在不影响计算精度的前提下有效地将这些非必要项剔除就显得尤为重要。基于逐步回归进行响应面基函数的显著性分析，其基本思想为通过多次引入、检验和剔除，以保证最终模型中的参数全部显著。显著性检验判定准则为

(5)

式(5)中：SSEh、SSEh+1分别为h、h+1项基函数RSM模型的响应误差平方和；κh、κh+1分别为h、h+1项基函数RSM模型的自由度。

对于给定的显著性水平α，当方差分析统计量Fh+1>Fα(1,N-p-2)时，则判定第h+1项显著，需要将该项引入到RSM模型中，其中p为响应面函数中非常数项项数。

1.3 响应面精度检验

(6)

(7)

(8)

2 基于SSA的结构有限元模型修正

2.1 麻雀搜索算法

麻雀搜索算法是模仿麻雀觅食行为和反哺行为而提出的一种新型群智能仿生优化算法，其大致流程如下。

(1)种群初始化。输入初始化麻雀种群数及相关系数。在n维全局变量空间中，m只麻雀构成的种群空间位置为

(9)

输出当前麻雀适应度值并进行顺序排列，将初始种群中适应度值较好的个体定义为发现者，较差的为跟随者。

(2)发现者位置更新。在觅食区域周围无被捕食风险时，发现者可开启广泛的搜索模式，引导种群向更高适应度值逼近；当种群边缘麻雀发现捕食者，迅速发出报警信号，当预警值大于安全值，种群立即做出反捕食行为，发现者将所有跟随者带领到安全区域觅食。发现者的位置更新为

(10)

式(10)中：i、j分别为麻雀数和参数维度，i=1,2,…,m，j=1,2,…,n；itermax为算法最大迭代次数；Xij为第i只麻雀在第j维中的位置信息；U为服从标准正态分布的随机数；T为单位行向量；δ为(0, 1]的随机数；WV(WV∈[0, 1])、SV(SV∈[0.5, 1])分别为预警值和安全值。

(3)跟随者位置更新。发现者和跟随者的身份是动态交互的。跟随者的能量越低，处在种群中的位置就越差，为获取更高能量，它们总能够找到能提供最好资源的发现者，通过围绕在发现者周围获取食物；跟随者也可通过竞争方式争夺食物。跟随者的位置更新为

(11)

式(11)中：XW、XP分别为当前种群中麻雀的最差位置和发现者最优位置；ω=rand{-1, 1}。

(4)侦察预警并更新位置。处于种群边缘的麻雀，是捕食者的首要目标，极易受到猎食。当边缘麻雀意识到危险时，会迅速向安全位置靠拢，进行相应的位置更新；位于群体中间位置的麻雀也会意识到危险迅速向其他麻雀随机靠近以躲避捕食者。位置更新公式为

(12)

式(12)中：XB为麻省的最优位置；β为步长控制参数，服从均值为0、标准差为1的正态分布随机数；ξ为[-1, 1]均匀随机数；ε为无限接近零的常数；fi为当前麻雀的适应度值；fg、fw分别为当前全局最优和最差适应度值。

2.2 有限元模型修正流程

基于RSM-SSA的结构有限元模型修正流程如图1所示。其基本流程可划分如下。

图1 基于RSM-SSA的有限元模型修正流程图

(1)明确目标函数，基于响应面法确定结构响应与参数之间近似的映射关系，生成有限元模型结果与实际值的残差和函数。

(2)运行MATLAB软件编制的SSA算法程序进行待修正参数的迭代优化。

(3)获取目标函数修正后的的全局最优解，修正后的有限元模型与实际响应再进行对比，以实际响应验证有限元修正结果。

3 数值算例

3.1 数值模型概况

以图2所示的Euler-Bernoulli梁为研究对象，基于SSA对其进行模型修正，验证算法的精度及效率。该梁长3 m，截面尺寸0.2 m×0.15 m，材料弹性模量为E=3.2×104MPa，材料密度取ρ=2 500 kg/m3。采用空间梁单元建立有限元模型，共划分10个单元，11个节点。在数值算例中，考虑对部分截面刚度折减的方式模拟实际梁体在运营过程中的局部损伤，以无损伤梁的有限元模型作为初始模型，局部损伤梁的有限元模型作为实际模型，并将局部损伤梁模型的模态自振频率作为梁动力特性的实测数据，从而达到模拟试验的结果。

图2 悬臂梁数值模型

假定梁体在单元2、3、5、6和8出现局部损伤，其截面刚度分别下降了初始值的10%、18%、24%、20%和15%。基于有限元计算分析，可以得到初始模型和实际模型的前5阶自振频率，如表1所示。两种模型前5阶自振频率误差均集中在5%附近，如果在实际工程结构中，初始有限元模型与实际值误差能够达到如此精度，可不进行模型修正，但对简单悬臂梁数值模型来说，不存在复杂的简化假设和建模误差，更不存在工程建设过程中各种不确定性因素带来的误差。显然，5.621%的最大误差对简单悬臂梁来说是无法接受的，所以为了满足有限元模型能较好地契合实际结构响应的需求，要对初始有限元模型进行修正，以获取精确的基准有限元模型。

表1 初始模型与实际模型自振频率对比

选择10个单元的刚度折减系数x(x=实际刚度/初始刚度)为待修正参数，即无损伤梁的初始参数向量为单位向量[1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0]T；局部损伤梁修正后的参数向量解为[1.0, 0.9, 0.82, 1.0, 0.76, 0.8, 1.0, 0.85, 1.0, 1.0]T。以结构前5阶自振频率为修正目标，利用SSA进行有限元模型的修正，修正后的有限元模型动力特性和实际模型响应应保持一致。

3.2 响应面构建

考虑到真实梁体在服役过程中由于材料老化、荷载作用及环境影响等，其使用性能呈现弱化趋势，故将各单元的刚度折减系数区间取为[0.7, 1.0]，并以材料弹性模量的折减反映刚度折减。基于LHD进行10因素30水平样本点抽取，得到30组样本数据，代入有限元模型中计算各试验水平的前5阶自振频率。部分LHD样本点与响应值如表2、表3所示。

表2 LHD样本点

表3 结构响应值

基于逐步回归理论，采用F检验分析待修正参数对响应函数的显著性。取显著性水平为0.05，由方差统计分析计算各参数项及其二次项的显著性水平F值，结果如图3、图4所示。

图3 各参数项对响应函数显著性分析

图4 各参数二次项对响应函数显著性分析

根据参数显著性分析结果，忽略显著性较低的参数项，在不影响精度的条件下达到优化方程的目的。采用不含交叉项的二次多项式对样本数据集进行响应面回归，实现结构模态自振频率与各单元刚度折减系数间复杂隐式关系的解析式表达，拟合后的响应面方程如式(13)～式(17)所示。响应面模型精度检验结果如表4所示。

(13)

(14)

(15)

(16)

(17)

由表4可知，5个响应面模型的决定系数均达到了0.99以上，说明该响应面方程与有限元计算结果吻合较好，拟合优度较高，客观真实地反映了结构频率与各参数之间的关系，具有足够的精度代替有限元模型进行计算分析。

表4 响应面模型精度检验

3.3 有限元模型修正

构建数学优化问题的目标函数为

Y(x1,x2,…,x10)=(Y1-9.058)2+(Y2-

12.064)2+(Y3-55.009)2+(Y4-

72.821)2+(Y5-153.064)2

(18)

根据式(1)将模型修正问题转化为有约束目标函数极值求解问题，得到如式(19)所示的数学模型为

(19)

利用MATLAB通用数学软件编制麻雀搜索算法运行程序，求解式(19)的数学模型。同时，为验证本文算法对目标函数的优化精度和效率，将其结果与其他新兴群智能优化算法如灰狼优化算法(grey wolf optimizer, GWO)、鲸鱼优化算法(whale optimization algorithm, WOA)进行对比。在SSA算法程序中：初始种群数目为20，预警值取0.8，发现者和侦查者比例取20%和10%，最大迭代次数100次，搜索过程中各变量的上下范围在[0.7, 1.0]，以10次优化结果的平均值作为全局最优解。为消除算法测试误差的影响，三种算法设置相同的初始种群数和迭代次数。基于SSA、GWO和WOA的模型修正结果如表5所示。自振频率修正后的误差与初始误差结果对比如图5所示。

由表5可知，对比基于不同算法修正后各单元刚度折减系数特征值与精确解结果显示，采用GWO、WOA修正得到的参数修正值与实际值最大相对误差分别为15.380%、15.714%，误差均值为9.043%、9.539%；而采用SSA修正计算的最大相对误差为12.463%，且误差均值仅为6.549%，参数的修正效果显著提升。修正的前5阶自振频率误差均值在三种算法中也是最小的，除悬臂梁5阶自振频率外，其余4阶自振频率误差均有明显降低。说明基于SSA的有限元模型修正方法比基于GWO、WOA的修正方法的精度和效率更高。从图5可以看出，修正前结构前5阶实际自振频率与初始有限元模型计算结果相差较大，最大误差为5.621%，其余各阶频率也有不小的误差。采用本文SSA修正得到的基准有限元模型响应值与实际值误差为0.183%～0.398%，表明修正后的模型能较好的反映结构真实响应，具有较高的精度。

图5 有限元模型初始误差与修正误差

表5 模型修正结果对比

4 结论

首次提出基于RSM-SSA的结构有限元模型修正基本方法，详尽阐述了该方法用于有限元模型修正的一般流程。以某高维局部损伤悬臂梁模型为数值算例，验证方法的可靠性和优越性。得到了以下结论。

(1)算法的合理选择是解决模型修正问题的基础。基于SSA的有限元模型修正方法获得了较好的修正结果，修正后模型各阶频率与实际频率最大误差仅为0.398%。修正精度和效率较其他群智能优化算法有显著提升，能明显降低模型响应与实际响应误差，真实地反映结构实际力学行为，修正后的基准有限元模型可为结构状态评估、系统识别和健康监测提供有效的分析手段。

(2)在处理高维问题时，待修正参数众多，结构响应与参数间往往呈现出强非线性，进行参数显著性分析可将灵敏度较弱的参数忽略以获得最优响应面，进一步提升修正效率。而且，可靠的试验设计方法将直接影响到响应面替代模型的拟合精度与试验成本。采用的拉丁超立方设计很好地处理了精度与试验成本这一矛盾问题，以较小的样本反映总体的变异规律，简单而不失精度。

(3)面对日趋复杂的实际工程结构，群智能优化算法作为人工智能领域的重要分支展现出极大的优势，未来其发展将从简单模拟低等生物种群向高等生物快速迈进，将在科学研究和工程建设领域都有很广阔的应用前景。