结构加速度频响函数模型修正的Kriging方法

工程科学学报，第39卷，第7期：1087-1093,2017年7月 Chinese Journal of Engineering,Vol.39,No.7:1087-1093,July 2017 D0l:10.13374/j.issn2095-9389.2017.07.015:http://journals..ustb.edu.cn 结构加速度频响函数模型修正的Kriging方法 王巨涛”，王春洁2)⑧，宋顺广》 1)北京航空航天大学机械工程及自动化学院，北京1001912)北京航空航天大学虚拟现实技术与系统国家重点实验室，北京100191 3)北京军立方机器人科技有限公司，北京100176 ☒通信作者，E-mail:wangcj(@buaa.edu.cn 摘要为提高结构频响函数模型修正效率，提出将Kigg模型引入优化过程，代替有限元模型进行迭代运算.基于频响曲 线对应频率点处的响应值之差构造目标函数，并结合初选设计参数进行实验设计.根据实验设计结果进行各参数的灵敏度 分析，进而筛选出模型修正的待修正参数，基于该参数及其响应构造Kriging模型，经检验有效的Kriging模型将参与模型修正 过程.以GARTEUR飞机模型为算例，基于加速度频响数据进行模型修正，修正后模型不仅能复现检验点处频响曲线，还能成 功预测结构局部修改后的频响曲线，证明了Kigg方法应用于频响函数模型修正的有效性. 关键词加速度频响函数；模型修正：Kriging;参数筛选 分类号V221·.92:TH113.1 Structural model updating of acceleration frequency response function based on Kriging method WANG Ju-tao》,WANG Chun-ie2a,SONG Shun--guang》 1)School of Mechanical Engineering and Automation,Beihang University,Beijing 100191,China 2)State Key Laboratory of Virtual Reality Technology and Systems,Beihang University,Beijing 100191,China 3)Beijing Junlifang Robot Technology Co.Ltd.Beijing 100176,China Corresponding author,E-mail:wangej@buaa.edu.cn ABSTRACT To improve the efficiency of frequency response function (FRF)model updating,the Kriging model is introduced into its optimization procedure,thus replacing the finite element model for iteration.The objective function is established based on the error of FRF curve at each corresponding frequency point,and a design of experiment (DOE)of the primary design parameters is per- formed.The sensitivity of each design is then assessed based on DOE results,and the parameters involved in model updating are se- lected.The Kriging model is constructed based on these selected parameters,and after confirming its response it participates in the model updating process.The GARTEUR aircraft model is used as an example,and it is found that based on acceleration FRF data, the FRF curve reappears at the checking point after model updating.The FRF curve is also predicted successfully when the structure is changed locally,and therefore the effectiveness of the Kriging method applied to FRF model updating is demonstrated. KEY WORDS acceleration frequency response function:model updating:Kriging:parameter selection 近年来，在结构动力学模型修正中，以基于模态参型建模误差还大0，而频响函数型修正方法避开了结 数和基于频响函数的修正方法最为常见，后者可直接构模态分析，适用于模态分布较密集的模型修正回，且 利用振动试验测得的频响数据进行修正而具有一定的 频响函数能够提供更多的数据，每条曲线都可以作为 便利性.在某些情况下，模态参数的识别比有限元模 目标函数来修正四.文献45]分别从理论上介绍了 收稿日期：201607-14 基金项目：国家自然科学基金重点资助项目(51635002)

工程科学学报，第 39 卷，第 7 期: 1087--1093，2017 年 7 月 Chinese Journal of Engineering，Vol． 39，No． 7: 1087--1093，July 2017 DOI: 10． 13374 /j． issn2095--9389． 2017． 07． 015; http: / /journals． ustb． edu． cn 结构加速度频响函数模型修正的 Kriging 方法 王巨涛1) ，王春洁1，2) ，宋顺广3) 1) 北京航空航天大学机械工程及自动化学院，北京 100191 2) 北京航空航天大学虚拟现实技术与系统国家重点实验室，北京 100191 3) 北京军立方机器人科技有限公司，北京 100176 通信作者，E-mail: wangcj@ buaa． edu． cn 摘 要 为提高结构频响函数模型修正效率，提出将 Kriging 模型引入优化过程，代替有限元模型进行迭代运算． 基于频响曲 线对应频率点处的响应值之差构造目标函数，并结合初选设计参数进行实验设计． 根据实验设计结果进行各参数的灵敏度 分析，进而筛选出模型修正的待修正参数，基于该参数及其响应构造 Kriging 模型，经检验有效的 Kriging 模型将参与模型修正 过程． 以 GAＲTEUＲ 飞机模型为算例，基于加速度频响数据进行模型修正，修正后模型不仅能复现检验点处频响曲线，还能成 功预测结构局部修改后的频响曲线，证明了 Kriging 方法应用于频响函数模型修正的有效性． 关键词 加速度频响函数; 模型修正; Kriging; 参数筛选 分类号 V221 + . 92; TH113. 1 Structural model updating of acceleration frequency response function based on Kriging method WANG Ju-tao1) ，WANG Chun-jie1，2)  ，SONG Shun-guang3) 1) School of Mechanical Engineering and Automation，Beihang University，Beijing 100191，China 2) State Key Laboratory of Virtual Ｒeality Technology and Systems，Beihang University，Beijing 100191，China 3) Beijing Junlifang Ｒobot Technology Co． Ltd． ，Beijing 100176，China  Corresponding author，E-mail: wangcj@ buaa． edu． cn ABSTＲACT To improve the efficiency of frequency response function ( FＲF) model updating，the Kriging model is introduced into its optimization procedure，thus replacing the finite element model for iteration． The objective function is established based on the error of FＲF curve at each corresponding frequency point，and a design of experiment ( DOE) of the primary design parameters is per￾formed． The sensitivity of each design is then assessed based on DOE results，and the parameters involved in model updating are se￾lected． The Kriging model is constructed based on these selected parameters，and after confirming its response it participates in the model updating process． The GAＲTEUＲ aircraft model is used as an example，and it is found that based on acceleration FＲF data， the FＲF curve reappears at the checking point after model updating． The FＲF curve is also predicted successfully when the structure is changed locally，and therefore the effectiveness of the Kriging method applied to FＲF model updating is demonstrated． KEY WOＲDS acceleration frequency response function; model updating; Kriging; parameter selection 收稿日期: 2016--07--14 基金项目: 国家自然科学基金重点资助项目( 51635002) 近年来，在结构动力学模型修正中，以基于模态参 数和基于频响函数的修正方法最为常见，后者可直接 利用振动试验测得的频响数据进行修正而具有一定的 便利性． 在某些情况下，模态参数的识别比有限元模 型建模误差还大［1］，而频响函数型修正方法避开了结 构模态分析，适用于模态分布较密集的模型修正［2］，且 频响函数能够提供更多的数据，每条曲线都可以作为 目标函数来修正［3］． 文献［4--5］分别从理论上介绍了

·1088· 工程科学学报，第39卷，第7期 设计参数型频响函数修正法和工程应用中的优势及局 A.(o）-A(o） 限.Esfandiari等采用基于频响函数的模型修正方 R,= (4) A(@:) 法，实现了对一标准桁架结构的损伤识别.文献了9] 式中，n,表示所选观测点数量，n表示在频响曲线上选 均采用加速度频响函数的模型修正方法，分别对某双 取的频率点个数，A.是试验模型在第i频率点的加速 星系统、桁架和简支梁的有限元模型进行模型修正 度幅值，A是有限元模型在第i频率点的加速度幅值， 工程结构有限元模型的动力学分析往往比较耗 w为第i个频率值. 时，这给基于优化的动力学模型修正的效率带来不小 2 Kriging法构造近似模型 困难，因此，能够实现快速迭代的近似模型方法越来越 受到关注，该方法通过建立结构输入（设计参数）与输 Kriging法构造的模型通过插值技术实现，属于线 出（结构响应）间的近似关系模型，用以取代有限元模 性回归分析的一种改进技术，被视为线性无偏估计 型参与优化过程，且利于采用不同的优化算法来研究 模型包含了线性回归部分和非参数部分： 问题，从而提高优化效率，这使该方法越来越具工程应 y(x)=F(B,x)+z(x)=f(x)B+z(x).(5) 用意义.刘洋等@分别介绍了两种常用的近似模型， 式中，B是回归系数，f(x)是样本点x的多项式函数， 并进行了基于模态参数的模型修正，指出了近似模型 z(x)为随机分布. 在结构有限元模型修正应用中的关键问题.文献1一 在设计空间中，f(x)作为模型构造中的全局近 12]均采用响应面法分别对基准模型模态和悬臂梁结 似，其形式类似于响应面法中的多项式.z(x)作为模 构多目标进行了模型修正，使修正后的各阶频率与试 型构造中局部偏差的近似，服从正态分布N(0,σ)， 验频率误差明显减小，文献3-l4]基于Kriging模型 其协方差不为零.随机分布z(x)的存在是Kriging模 分别对多体动力学模型和拱桥模型进行了修正· 型与传统响应面的主要不同点，z(x)的协方差矩 阵为： 加速度频响函数修正方法 Cov (x),z(x)]=2R. (6) 加速度频响函数是频响函数的一种，假设n自由 式中，σ，为z(x)的方差：x和x是所有训练样本点中的 度阻尼系统中，基本动力学方程为： 任意两个；R是由R,(x,x)组成的对称矩阵 Mi+Ci+Kx=F (1) R,是相关函数，该函数能够表征样本点之间的空 式中，M为系统质量阵，C为系统阻尼阵，K为系统刚 间相关性，直接影响模型的近似精度.R,的形式如下： 度阵，x为位移向量，F为系统激励 ,,)=ΠR(a,d),d,=ld-1.() 在简谐激励下，经傅里叶变换可得频域内输入和 输出的关系为： 式中，N,是设计变量的数量，0，为相关系数，x和x分 A(0)=H(@)F(o) (2) 别是x:和x样本点第k个分量. 式中，A(o)为稳态响应，H(o)为加速度频响函数， 本文的R,选取应用较广的高斯相关函数回，其核 F(a)为简谐激励. 函数的形式为： 可得有限元模型的加速度频响函数矩阵为： R(0,d)=exp(-0d) (8) H()=[M-ic 因σ和R都是，的函数，且o和R可通过响应的 w2」· (3) 最大化似然估计来确定，如式(9)所示，求解σ，和R 式中，ω为激励频率. 即为确定0值 在建立的有限元模型中，所用结构参数与试验模 N.In (o2)+In IRI 型结构参数间通常会存在误差，这直接或间接体现在 用 (9) 2 质量误差、阻尼误差和刚度误差中，从而使有限元模型 式中，N表示训练样本点个数 的频响函数不同于试验模型的响应.频响函数模型修 由于式(5)中的随机分布z(x)已经有足够的能力 正过程，即为利用试验测得的频响数据，通过调整有限 描述函数响应的变化趋势，其回归部分(x)可以视 元模型的部分结构参数（待修正参数）来减小与试验 为常数项，则0为Kriging模型的唯一未知数. 模型间响应误差 当Kriging模型y(x)构造完成，在待测点x。的响 模型修正过程中，除了确定有限元模型的待修正 应值(x。)为： 参数作为优化过程的设计变量外，还需构造优化过程 y(xo)=B'+r(xo)R-(Y-fB"). (10) 的目标函数.本文利用频响函数曲线在同一频率点处 其中： 加速度幅值的差值，构造出了单自由度下的目标函数 B=(fR-f)fR-Y. (11) 形式为： 式中，Y为训练样本点响应组成的列向量：r(xo)为训

工程科学学报，第 39 卷，第 7 期 设计参数型频响函数修正法和工程应用中的优势及局 限． Esfandiari 等［6］采用基于频响函数的模型修正方 法，实现了对一标准桁架结构的损伤识别． 文献［7--9］ 均采用加速度频响函数的模型修正方法，分别对某双 星系统、桁架和简支梁的有限元模型进行模型修正． 工程结构有限元模型的动力学分析往往比较耗 时，这给基于优化的动力学模型修正的效率带来不小 困难，因此，能够实现快速迭代的近似模型方法越来越 受到关注，该方法通过建立结构输入( 设计参数) 与输 出( 结构响应) 间的近似关系模型，用以取代有限元模 型参与优化过程，且利于采用不同的优化算法来研究 问题，从而提高优化效率，这使该方法越来越具工程应 用意义． 刘洋等［10］分别介绍了两种常用的近似模型， 并进行了基于模态参数的模型修正，指出了近似模型 在结构有限元模型修正应用中的关键问题． 文献［11-- 12］均采用响应面法分别对基准模型模态和悬臂梁结 构多目标进行了模型修正，使修正后的各阶频率与试 验频率误差明显减小，文献［13--14］基于 Kriging 模型 分别对多体动力学模型和拱桥模型进行了修正． 1 加速度频响函数修正方法 加速度频响函数是频响函数的一种，假设 n 自由 度阻尼系统中，基本动力学方程为: M x·· + C x· + Kx = F． ( 1) 式中，M 为系统质量阵，C 为系统阻尼阵，K 为系统刚 度阵，x 为位移向量，F 为系统激励． 在简谐激励下，经傅里叶变换可得频域内输入和 输出的关系为: A( ω) = H( ω) F( ω) ． ( 2) 式中，A( ω) 为稳态响应，H( ω) 为加速度频响函数， F( ω) 为简谐激励． 可得有限元模型的加速度频响函数矩阵为: H( ω) [ = M － jC ω － K ω2 ] － 1 ． ( 3) 式中，ω 为激励频率． 在建立的有限元模型中，所用结构参数与试验模 型结构参数间通常会存在误差，这直接或间接体现在 质量误差、阻尼误差和刚度误差中，从而使有限元模型 的频响函数不同于试验模型的响应． 频响函数模型修 正过程，即为利用试验测得的频响数据，通过调整有限 元模型的部分结构参数( 待修正参数) 来减小与试验 模型间响应误差． 模型修正过程中，除了确定有限元模型的待修正 参数作为优化过程的设计变量外，还需构造优化过程 的目标函数． 本文利用频响函数曲线在同一频率点处 加速度幅值的差值，构造出了单自由度下的目标函数 形式为: Ｒnp = ∑ np n = 1 ∑ nf i = 1 Ati ( ωi ) － Aai ( ωi ) Ati ( ωi ) ． ( 4) 式中，np表示所选观测点数量，nf表示在频响曲线上选 取的频率点个数，Ati是试验模型在第 i 频率点的加速 度幅值，Aai是有限元模型在第 i 频率点的加速度幅值， ωi为第 i 个频率值． 2 Kriging 法构造近似模型 Kriging 法构造的模型通过插值技术实现，属于线 性回归分析的一种改进技术，被视为线性无偏估计． 模型包含了线性回归部分和非参数部分: y( x) = F( β，x) + z( x) = f T ( x) β + z( x) ． ( 5) 式中，β 是回归系数，f( x) 是样本点 x 的多项式函数， z( x) 为随机分布． 在设计空间中，f( x) 作为模型构造中的全局近 似，其形式类似于响应面法中的多项式． z( x) 作为模 型构造中局部偏差的近似，服从正态分布 N( 0，σ2 z ) ， 其协方差不为零． 随机分布 z( x) 的存在是 Kriging 模 型与传 统 响 应 面 的 主 要 不 同 点，z ( x) 的 协 方 差 矩 阵为: Cov［z( xi ) ，z( xj ) ］= σ2 z Ｒ． ( 6) 式中，σz为 z( x) 的方差; xi和 xj是所有训练样本点中的 任意两个; Ｒ 是由 Ｒij( xi，xj ) 组成的对称矩阵． Ｒij是相关函数，该函数能够表征样本点之间的空 间相关性，直接影响模型的近似精度． Ｒij的形式如下: Ｒij( xi，xj ) = ∏ Nv k = 1 Ｒk ( θk，dk ) ，dk = | xk i － xk j | ． ( 7) 式中，Nv是设计变量的数量，θk为相关系数，xk i 和 xk j 分 别是 xi和 xj样本点第 k 个分量． 本文的 Ｒij选取应用较广的高斯相关函数［15］，其核 函数的形式为: Ｒk ( θk，dk ) = exp ( － θkd2 k ) ． ( 8) 因 σz和 Ｒ 都是 θk的函数，且 σz和 Ｒ 可通过响应的 最大化似然估计来确定［15］，如式( 9) 所示，求解 σz和 Ｒ 即为确定 θk值． max θk ＞ 0 ( － Nsln ( σ2 z ) + ln | Ｒ| ) 2 ． ( 9) 式中，Ns表示训练样本点个数． 由于式( 5) 中的随机分布 z( x) 已经有足够的能力 描述函数响应的变化趋势，其回归部分 f T ( x) 可以视 为常数项，则 θk为 Kriging 模型的唯一未知数． 当 Kriging 模型 y( x) 构造完成，在待测点 x0 的响 应值 y^( x0 ) 为: y^( x0 ) = β* + r T ( x0 ) Ｒ － 1 ( Y － fβ* ) ． ( 10) 其中: β* = ( f T Ｒ － 1 f) － 1 f T Ｒ － 1Y． ( 11) 式中，Y 为训练样本点响应组成的列向量; r( x0 ) 为训 · 8801 ·

王巨涛等：结构加速度频响函数模型修正的Kriging方法 ·1089* 练样本点与待测位置x。所组成的相关向量，其形式如 建模，整体采用壳单元，在模型中尾翼与机身的连接 下式所示 处，以及机翼与机身的连接处都设立了弹簧单元，作为 r(xo)=Ro(o,x1）,R2(xo,x2）,…,Rav(xo,xx)]T 刚度类设计参数.此外，考虑到零件间连接时所用螺 (12) 栓的附加质量，在有限元模型的相应位置增加了质 量单元，最终的有限元模型如图2所示，共包括2332 3模型修正流程 个节点，2010个壳单元，2个弹簧单元，4个质量点 根据前两节理论基础，提出将Kriging近似模型引 单元 入结构有限元模型修正过程，用该模型代替有限元模 图2中1和2所示为弹簧单元和质量单元所处位 型参与优化求解，达到频响函数模型修正的目的 置，所示点3为频响函数观测点，4为单位简谐激励位 构造Kriging模型之前先对初选参数进行实验设 置.对该结构进行2%模态阻尼下的无约束的自由频 计，基于实验设计结果分析各参数对输出响应的灵敏 响分析，选取点3Y自由度频响函数的加速度幅值之差 度，挑选出灵敏度较高的设计参数作为构造Kriging模 来构造目标函数，并令目标函数值作为本算例的输出 型的输入参数，再结合实验设计的输出响应结果构造 响应.为检验模型修正效果，设置所示点5为修正后 Kriging模型. 模型检验点.此外，为检验修正后模型的预测能力，所 对构造出的Kriging模型进行近似精度评估，当其 示I和Ⅱ为结构的局部质量修改位置 精度满足要求时可代替有限元模型进行优化求解，优 化目标为有限元模型与试验模型的响应误差最小.最 后利用优化所得设计参数的最优解更新有限元模型， 即可得修正后有限元模型，通过提取频响函数曲线查 看模型修正效果.本文提出的模型修正流程如图1 所示. 有限元 试险 模型 频响数据 初选设计参数 构造目标函数 试验设计 优化求解 图2 GARTEUR飞机结构有限元模型图 Fig.2 Finite element of GARTEUR aircraft structure 参数筛选 <满足误差要求> 4.1实验设计与参数筛选 本文采用图2所示初始有限元模型作为本算例的 构造Kriging模型 是 “试验模型”，相应的“有限元模型”通过对试验模型不 参数修正值 同类型的11个可能存在误差的设计参数值进行偏离 满足近似 精度要求 后获得，具体参数及偏离情况如表1所示，机翼受阻尼 修正后 有限元模型 材料的影响其材料属性单独作为一个设计参数，表中 5个参数值增加10%，其余6个参数值减小10%，各参 图1模型修正流程图 Fig.I Flowchart of model updating 数即为有限元模型的初始设计参数，该参数将作为实 验设计的因子. 4 算例验证 实验设计的响应对应频响函数幅值误差，对试验 模型和有限元模型分别进行2~100Hz内的自由频响 GARTEUR SM-AG-19飞机模型，是欧洲航空科 分析，载荷步为1Hz.可获得观测点的加速度频响函 技研究组织中结构与材料工作组建立的一个典型的标 数的实部曲线，根据该曲线的特征去掉峰值附近4~8 准模型.该飞机模型具有真实飞机的高柔度、模态频 Hz和65~73Hz的频率点，令其不参与目标函数值的 率低和模态密集的特点.模型主体由铝制结构组成， 计算叨，按式(4)构造的此次优化过程的目标函数如 特殊之处在于机翼上表面为含约束层的黏弹性阻尼 式(13)所示.考虑到实际数据测量中试验数据会受到 材料n▣ 噪声污染，则给试验模型的频响曲线加入了10%的随 本文使用软件Patran对该飞机结构进行了有限元 机噪声

王巨涛等: 结构加速度频响函数模型修正的 Kriging 方法 练样本点与待测位置 x0所组成的相关向量，其形式如 下式所示． r( x0 ) =［Ｒ01 ( x0，x1 ) ，Ｒ02 ( x0，x2 ) ，…，Ｒ0Ns ( x0，xNs ) ］T ． ( 12) 3 模型修正流程 根据前两节理论基础，提出将 Kriging 近似模型引 入结构有限元模型修正过程，用该模型代替有限元模 型参与优化求解，达到频响函数模型修正的目的． 构造 Kriging 模型之前先对初选参数进行实验设 计，基于实验设计结果分析各参数对输出响应的灵敏 度，挑选出灵敏度较高的设计参数作为构造 Kriging 模 型的输入参数，再结合实验设计的输出响应结果构造 Kriging 模型． 对构造出的 Kriging 模型进行近似精度评估，当其 精度满足要求时可代替有限元模型进行优化求解，优 化目标为有限元模型与试验模型的响应误差最小． 最 后利用优化所得设计参数的最优解更新有限元模型， 即可得修正后有限元模型，通过提取频响函数曲线查 看模型修正效果． 本文提出的模型修正流程如图 1 所示． 图 1 模型修正流程图 Fig． 1 Flowchart of model updating 4 算例验证 GAＲTEUＲ SM--AG--19 飞机模型，是欧洲航空科 技研究组织中结构与材料工作组建立的一个典型的标 准模型． 该飞机模型具有真实飞机的高柔度、模态频 率低和模态密集的特点． 模型主体由铝制结构组成， 特殊之处在于机翼上表面为含约束层的黏弹性阻尼 材料［16］． 本文使用软件 Patran 对该飞机结构进行了有限元 建模，整体采用壳单元，在模型中尾翼与机身的连接 处，以及机翼与机身的连接处都设立了弹簧单元，作为 刚度类设计参数． 此外，考虑到零件间连接时所用螺 栓的附加质量，在有限元模型的相应位置增加了质 量单元，最终的有限元模型如图 2 所示，共包括 2332 个节点，2010 个 壳 单 元，2 个 弹 簧 单 元，4 个 质 量 点 单元． 图 2 中 1 和 2 所示为弹簧单元和质量单元所处位 置，所示点 3 为频响函数观测点，4 为单位简谐激励位 置． 对该结构进行 2% 模态阻尼下的无约束的自由频 响分析，选取点 3Y 自由度频响函数的加速度幅值之差 来构造目标函数，并令目标函数值作为本算例的输出 响应． 为检验模型修正效果，设置所示点 5 为修正后 模型检验点． 此外，为检验修正后模型的预测能力，所 示Ⅰ和Ⅱ为结构的局部质量修改位置． 图 2 GAＲTEUＲ 飞机结构有限元模型图 Fig． 2 Finite element of GAＲTEUＲ aircraft structure 4. 1 实验设计与参数筛选 本文采用图 2 所示初始有限元模型作为本算例的 “试验模型”，相应的“有限元模型”通过对试验模型不 同类型的 11 个可能存在误差的设计参数值进行偏离 后获得，具体参数及偏离情况如表 1 所示，机翼受阻尼 材料的影响其材料属性单独作为一个设计参数，表中 5 个参数值增加 10% ，其余 6 个参数值减小 10% ，各参 数即为有限元模型的初始设计参数，该参数将作为实 验设计的因子． 实验设计的响应对应频响函数幅值误差，对试验 模型和有限元模型分别进行 2 ～ 100 Hz 内的自由频响 分析，载荷步为 1 Hz． 可获得观测点的加速度频响函 数的实部曲线，根据该曲线的特征去掉峰值附近 4 ～ 8 Hz 和 65 ～ 73 Hz 的频率点，令其不参与目标函数值的 计算［17］，按式( 4) 构造的此次优化过程的目标函数如 式( 13) 所示． 考虑到实际数据测量中试验数据会受到 噪声污染，则给试验模型的频响曲线加入了 10% 的随 机噪声． · 9801 ·

·1090. 工程科学学报，第39卷，第7期 表1设计参数初始误差设置 型的训练样本x的数据来源，每个训练样本点只包含 Table 1 Initial error of design parameters 筛选出的3个参数的水平值及其对应的响应值，进而 参数 名称 初始误差/% 求解Kriging模型的构造参数. 尾翼x向弹簧刚度 -10 本文选取高斯相关函数作为训练样本点的空间相 42 尾翼y向弹簧刚度 0 关函数，通过求解式(9)可得修正参数“材料密度”、 A3 尾翼Σ向弹簧刚度 -10 “弹性模量”和“机翼模量”所对应的日值分别为：6，= A 机翼x向弹簧刚度 10 4.95259,02=4.45985,03=5.33281.则相应的B、R As 机翼y向弹簧刚度 -10 和r随之确定，Kriging模型便构造完成，可以预测设计 6 机翼：向弹簧刚度 10 空间中任一点的输出响应.当其近似精度满足要求便 M 材料密度 -10 认为模型有效，即可代替有限元模型进行优化求解. As 弹性模量 10 为评估近似模型的有效性，用实验设计的方法选 取训练样本点之外的l0组设计参数值，作为Kriging Ag 机翼质量单元质量 -10 Ao 模型近似精度评估点.通过计算该10个评估点处有 尾翼质量单元质量 10 限元计算值和Kriging模型预测值间的决定系数 Au 机翼模量 -10 (COD)和均方根误差(RMSE)来定量分析近似精度， 注：机翼模量表示机翼材料的弹性模量 当COD越接近l或RMsE越接近0时表明Kriging模 |A.（w）-A（o） 型预测越精确.COD和RMSE分别通过式(14)和式 (13) A.(o） (15)计算，所选的10个评估点处的响应值如图4 在各因子初始值的20%区间内进行实验设计，用 所示. 采样较为均匀的优化拉丁方设计生成因子的设计矩 阵，并最终产生了200个样本点.根据该样本点的因 (y-)2 COD=1- 女 (14) 子和响应数据建立多元二次回归方程，从而计算各因 ∑ (y-2 子对响应的灵敏度值，该值经过归一化并转换为百分 比的形式，按从大到小排列如图3所示 RMSE=√ (y-)2. (15) 60 式中，k为评估点数量，y为评估点有限元模型响应的 50 均值，y,为评估点有限元模型响应值，：为Kriging模型 在评估点处的预测值. 30 250m 。有限元 230 -◆-Kriging 20 210 190 10 170 0 ☑☑man山 150 AA,AA:小AAA,AA。Am 130 因子 110 图3各因子灵敏度 90 Fig.3 Sensitivity of each factor 70 50 由图3可见，因子A(机翼模量)和A,(材料密 5 6 7 8 910 评估点 度)的灵敏度最高，在优化过程中肯定能得到有效优 化.此外，经简单尝试把灵敏度第三高的因子A,(弹性 图4各评估点处响应值 Fig.4 Response value of each assessment point 模量)与前两个因子一同优化时，因子A也能得到可 接受的优化效果，而更低灵敏度因子便不能有效优化， 经计算，此次修正所用Kriging模型COD值和 甚至出现优化后初始误差被扩大的现象.因此，选择 RMSE值分别约为0.984和0.0216，表明其对有限元 因子A,、A,和A作为待修正参数参与优化过程 模型具有良好的近似精度，该模型有效 4.2 Kriging模型的构造及评估 4.3模型修正 基于筛选出的3个参数来构造Kriging近似模型， 本文采用模拟退火法对模型修正的优化问题进行 需要将实验设计生成的200个样本点作为Kriging模 求解，该算法类似于遗传算法，但比遗传算法的编制简

工程科学学报，第 39 卷，第 7 期 表 1 设计参数初始误差设置 Table 1 Initial error of design parameters 参数 名称 初始误差/% A1 尾翼 x 向弹簧刚度 － 10 A2 尾翼 y 向弹簧刚度 10 A3 尾翼 z 向弹簧刚度 － 10 A4 机翼 x 向弹簧刚度 10 A5 机翼 y 向弹簧刚度 － 10 A6 机翼 z 向弹簧刚度 10 A7 材料密度 － 10 A8 弹性模量 10 A9 机翼质量单元质量 － 10 A10 尾翼质量单元质量 10 A11 机翼模量 － 10 注: 机翼模量表示机翼材料的弹性模量． Ｒ = ∑ 85 i = 1 Ati ( ωi ) － Aai ( ωi ) Ati ( ωi ) ． ( 13) 在各因子初始值的!20% 区间内进行实验设计，用 采样较为均匀的优化拉丁方设计生成因子的设计矩 阵，并最终产生了 200 个样本点． 根据该样本点的因 子和响应数据建立多元二次回归方程，从而计算各因 子对响应的灵敏度值，该值经过归一化并转换为百分 比的形式，按从大到小排列如图 3 所示． 图 3 各因子灵敏度 Fig． 3 Sensitivity of each factor 由图 3 可见，因子 A11 ( 机翼模量) 和 A7 ( 材料密 度) 的灵敏度最高，在优化过程中肯定能得到有效优 化． 此外，经简单尝试把灵敏度第三高的因子 A8 ( 弹性 模量) 与前两个因子一同优化时，因子 A8 也能得到可 接受的优化效果，而更低灵敏度因子便不能有效优化， 甚至出现优化后初始误差被扩大的现象． 因此，选择 因子 A7、A8和 A11作为待修正参数参与优化过程． 4. 2 Kriging 模型的构造及评估 基于筛选出的 3 个参数来构造 Kriging 近似模型， 需要将实验设计生成的 200 个样本点作为 Kriging 模 型的训练样本 xi的数据来源，每个训练样本点只包含 筛选出的 3 个参数的水平值及其对应的响应值，进而 求解 Kriging 模型的构造参数． 本文选取高斯相关函数作为训练样本点的空间相 关函数，通过求解式( 9) 可得修正参数“材料密度”、 “弹性模量”和“机翼模量”所对应的 θ 值分别为: θ1 = 4. 95259，θ2 = 4. 45985，θ3 = 5. 33281． 则相应的 β* 、Ｒ 和 r 随之确定，Kriging 模型便构造完成，可以预测设计 空间中任一点的输出响应． 当其近似精度满足要求便 认为模型有效，即可代替有限元模型进行优化求解． 为评估近似模型的有效性，用实验设计的方法选 取训练样本点之外的 10 组设计参数值，作为 Kriging 模型近似精度评估点． 通过计算该 10 个评估点处有 限 元 计 算 值 和 Kriging 模型预测值间的决定系数 ( COD) 和均方根误差( ＲMSE) 来定量分析近似精度， 当 COD 越接近 1 或 ＲMSE 越接近 0 时表明 Kriging 模 型预测越精确． COD 和 ＲMSE 分别通过式( 14) 和式 ( 15) 计 算，所 选 的 10 个 评 估 点 处 的 响 应 值 如 图 4 所示． COD = 1 － ∑ k i = 1 ( yi － y^ i ) 2 ∑ k i = 1 ( yi － y) 2 ． ( 14) ＲMSE = 1 k y ∑ k i = 1 ( yi － y^ i ) 槡 2 ． ( 15) 式中，k 为评估点数量，y 为评估点有限元模型响应的 均值，yi为评估点有限元模型响应值，y^ i 为 Kriging 模型 在评估点处的预测值． 图 4 各评估点处响应值 Fig． 4 Ｒesponse value of each assessment point 经计算，此 次 修 正 所 用 Kriging 模 型 COD 值 和 ＲMSE 值分别约为 0. 984 和 0. 0216，表明其对有限元 模型具有良好的近似精度，该模型有效． 4. 3 模型修正 本文采用模拟退火法对模型修正的优化问题进行 求解，该算法类似于遗传算法，但比遗传算法的编制简 · 0901 ·

王巨涛等：结构加速度频响函数模型修正的Kriging方法 ·1091· 单，结合概率突跳性使目标函数不至于陷入局部最优 不理想 解，能有效搜索全局最优解 根据表2中Kriging法的修正结果，对有限元模型 优化过程的设计变量对应3个待修正设计参数， 的相应设计参数值进行更新，保留其余8个未参与修 设计变量的取值范围为其初始值的20%.令4.2节构 正参数的初始值，则可得修正前后观测点3的Y自由 造的Kriging模型代替有限元模型参与优化迭代，经过 度频响函数的加速度幅值（实部）曲线，如图5所示， 最大1000次迭代后，得出使目标函数值最小的3个参 修正前约7Hz处的峰值误差，以及60Hz处的频率误 数的最优解即为修正值，如表2所示.至此，模型修正 差，这两处主要误差均得到了较好修正，证明了Kig 的求解过程已完成，总用时为200样本点生成时间和 ing法的有效性. Kriging模型的l000次迭代，相对于前者，后者的用时 6000 修正后 可以忽略不计.若采用该模拟退火算法直接迭代有限 5000 一试验值 元模型求解参数最优值，其运算时间将是本文所提修 4000 ·修正前 分 3000 正方法的5倍. 2000 此外，为了与响应面法构造近似模型时优化效果 1000 进行对比，基于相同的样本点构造了比较常用的2阶 响应面以及更高阶的3阶响应面，采用该两个响应面 -1000 修正后的参数值及误差如表2所示 -2000 -3000 表2筛选出参数的修正结果 400 Table 2 Updated result of selected parameters 0 102030405060708090100 频率Hz 材料密度/弹性模量/机翼模量/ 名称 (kg'm-2) GPa GPa 图5加速度频响曲线修正前后对比 修正前 Fig.5 Comparison of acceleration FRF curves 2430 78.1 63.9 试验值 2700 之 71 4.4 修正后模型的检验 Kriging法修正值 2676.2 68.234 69.915 为进一步检验Kriging法的修正效果，对观测点3 2阶响应面法修正值 2889.9 68.496 76.680 的Y自由度频响函数的虚部曲线进行修正前后对比， 3阶响应面法修正值 2916.0 74.218 76.680 如图6(a)所示.同时，对检验点5的Z自由度频响函 Kriging法修正后误差/% 0.88 3.9 1.5 数的实部曲线进行了修正前后对比，如图6(b)所示. 2阶响应面法修正后误差/% 7.0 3.5 8.0 可见，两个检验点处修正后有限元模型频响曲线与试 3阶响应面法修正后误差/% 8.0 4.5 8.0 验频响曲线基本重合 最后，为了检验修正后有限元模型的预测能力，对 由表2可见，当近似模型采用响应面法构造时，3 试验模型做如下修改：图2中I位置为机尾局部，在平 个设计参数的误差，除了2阶响应面模型的“弹性模 尾上增加一个0.92kg的集中质量：Ⅱ位置为机翼局 量”优化后误差稍小于Kriging法外，其余修正效果均 部，将原先翼尖上0.15kg的质量块替换为0.72kg.有 1000 60 (b) …修正后 一试验值 修正前 -1000 修正后 -2000 一试验值 -一修正前 -3000 4000 -5000 -40 -6000 60 -7000 102030405060 80 70 80 90 0 102030405060708090100 频率Hz 频率/Hz 图6结构变化后测点3和检验点5加速度频响曲线.()测点3加速度频响曲线虚部修正前后对比：(b)检验点5加速度颊响曲线修正 前后对比 Fig.6 Comparison of acceleration FRF curves at measuring Point 3 and checking Point 5:(a)comparison of imaginary parts before and after the up- dating at measuring Point 3:(b)comparison of real parts before and after the updating at checking Point 5

王巨涛等: 结构加速度频响函数模型修正的 Kriging 方法 单，结合概率突跳性使目标函数不至于陷入局部最优 解，能有效搜索全局最优解． 优化过程的设计变量对应 3 个待修正设计参数， 设计变量的取值范围为其初始值的!20% ． 令 4. 2 节构 造的 Kriging 模型代替有限元模型参与优化迭代，经过 最大 1000 次迭代后，得出使目标函数值最小的 3 个参 数的最优解即为修正值，如表 2 所示． 至此，模型修正 的求解过程已完成，总用时为 200 样本点生成时间和 Kriging 模型的 1000 次迭代，相对于前者，后者的用时 可以忽略不计． 若采用该模拟退火算法直接迭代有限 元模型求解参数最优值，其运算时间将是本文所提修 正方法的 5 倍． 此外，为了与响应面法构造近似模型时优化效果 进行对比，基于相同的样本点构造了比较常用的 2 阶 响应面以及更高阶的 3 阶响应面，采用该两个响应面 修正后的参数值及误差如表 2 所示． 表 2 筛选出参数的修正结果 Table 2 Updated result of selected parameters 名称 材料密度/ ( kg·m － 3 ) 弹性模量/ GPa 机翼模量/ GPa 修正前 2430 78. 1 63. 9 试验值 2700 71 71 Kriging 法修正值 2676. 2 68. 234 69. 915 2 阶响应面法修正值 2889. 9 68. 496 76. 680 3 阶响应面法修正值 2916. 0 74. 218 76. 680 Kriging 法修正后误差/% 0. 88 3. 9 1. 5 2 阶响应面法修正后误差/% 7. 0 3. 5 8. 0 3 阶响应面法修正后误差/% 8. 0 4. 5 8. 0 图 6 结构变化后测点 3 和检验点 5 加速度频响曲线． ( a) 测点 3 加速度频响曲线虚部修正前后对比; ( b) 检验点 5 加速度频响曲线修正 前后对比 Fig． 6 Comparison of acceleration FＲF curves at measuring Point 3 and checking Point 5: ( a) comparison of imaginary parts before and after the up￾dating at measuring Point 3; ( b) comparison of real parts before and after the updating at checking Point 5 由表 2 可见，当近似模型采用响应面法构造时，3 个设计参数的误差，除了 2 阶响应面模型的“弹性模 量”优化后误差稍小于 Kriging 法外，其余修正效果均 不理想． 根据表 2 中 Kriging 法的修正结果，对有限元模型 的相应设计参数值进行更新，保留其余 8 个未参与修 正参数的初始值，则可得修正前后观测点 3 的 Y 自由 度频响函数的加速度幅值( 实部) 曲线，如图 5 所示， 修正前约 7 Hz 处的峰值误差，以及 60 Hz 处的频率误 差，这两处主要误差均得到了较好修正，证明了 Krig￾ing 法的有效性． 图 5 加速度频响曲线修正前后对比 Fig． 5 Comparison of acceleration FＲF curves 4. 4 修正后模型的检验 为进一步检验 Kriging 法的修正效果，对观测点 3 的 Y 自由度频响函数的虚部曲线进行修正前后对比， 如图 6( a) 所示． 同时，对检验点 5 的 Z 自由度频响函 数的实部曲线进行了修正前后对比，如图 6( b) 所示． 可见，两个检验点处修正后有限元模型频响曲线与试 验频响曲线基本重合． 最后，为了检验修正后有限元模型的预测能力，对 试验模型做如下修改: 图 2 中Ⅰ位置为机尾局部，在平 尾上增加一个 0. 92 kg 的集中质量; Ⅱ位置为机翼局 部，将原先翼尖上 0. 15 kg 的质量块替换为 0. 72 kg． 有 · 1901 ·