采用两步式模型修正过程识别海上风电支撑结构的损伤

2016-10-12李英超王树青

海洋工程 2016年3期

关键词：正则振型修正

李英超，张敏，王树青

(1.鲁东大学土木工程学院，山东烟台 264025；2.中国海洋大学工程学院，山东青岛 266100)

采用两步式模型修正过程识别海上风电支撑结构的损伤

李英超1，张敏2，王树青2

(1.鲁东大学土木工程学院，山东烟台 264025；2.中国海洋大学工程学院，山东青岛 266100)

针对海上风电支撑结构损伤识别的难点，将一种新的模型修正方法应用到了三脚架式海上风电支撑结构的损伤识别中。采用数学上常用的正则化方法，解决了实测模态含噪声条件下损伤识别“病态”方程组的求解问题；结合Guyan扩阶技术，实现了基于低阶不完备模态的损伤识别。另外，为提高损伤识别精度，结合算例提出一种两步式模型修正过程：第一步，初步定位损伤位置；第二步，精确识别损伤位置和损伤程度。数值研究结果表明，采用两步式模型修正过程的方法对于海上风电支撑结构的损伤识别是有效的。

损伤识别；模型修正；海上风电；正则化；模态扩阶

Abstract:Through a numerical example,a recently developed model updating method is applied to identify the damages of the tripod-type offshore wind turbine support structure.Regularization method of mathematics is used to solve the “ill-conditioned” damage identification equations when noise polluted measured modes are used.With Guyan expansion technology,damage identification is realized based on a few low order incomplete modes.Moreover,in order to improve the accuracy of damage identification,a two-step model updating process is proposed.The first step is to locate damage position.The second step is to identify damage location and damage severity accurately.Results of numerical examples show that the method used in this paper is effective for damage identification of offshore wind structure.

Keywords:damage identification; model updating; offshore wind turbine; regularization; modal expansion

海上风电支撑结构长期作用在复杂的海洋环境中，除了要承受风机荷载外，还将受到波浪、流、风，甚至冰击、船的撞击和台风、地震等极端荷载的作用，另外结构本身还要遭到环境腐蚀、海生物附着、海底冲刷、复杂海床变化等影响。在此环境下，结构容易产生各种各样的损伤，使得结构整体动力特性发生改变，严重时还会导致结构失效，造成巨大的经济损失和不良的社会影响。为了合理有效地维护海上风电结构，减少因结构损伤导致的风机停机时间，降低维护成本，采用有效的方法对风电结构在服役期内的安全性进行实时监测或定期检测是非常必要的[1]。

然而，海上风电支撑结构大部分位于水面以下和泥面以下，且常以建筑群的形式存在，测试非常困难，传统的无损检测方法很难进行应用。与此相反，基于振动测试的结构整体检测技术日益受到国内外学者的关注，并成为研究的热点[2]。其基本原理是根据结构的实测动力特性反推结构的物理参数，从而诊断结构的健康性。其最突出的优点是利用环境激励下的动力响应进行损伤诊断，检测过程不影响正常使用[3]。

基于振动测试的结构损伤诊断技术在航空航天、土木工程和机械工程等领域得到了深入的研究和广泛的应用，目前学者们已经提出了多种较成熟的方法，如模型修正法、敏感性分析方法、神经网络方法、小波包变换法、遗传算法、频响函数法和模态应变能法等。但这些方法在海洋结构物中的应用还受到一些限制，如传感器布置困难，使得实测模态极度不完备；另外环境激励下的实测信号噪声污染较大，识别精度较差等。针对这类问题， Hu和Li等[4-5]提出了交叉模型交叉模态(cross model cross mode,简称CMCM)有限元模型修正方法，并将其应用于海洋平台结构的损伤诊断中。相比于传统方法，CMCM方法最大的优点是可以应用有限的低阶实测模态构建更多的线性方程，使得修正因子(损伤程度)更利于求解。然而这种方法在实际应用中存在两个技术瓶颈：1)实测模态空间不完备问题；2)实测模态含噪声条件下模型修正“病态”方程组的求解问题[6]。

针对第1)类问题，Li等[5]展开了研究，结合扩阶/缩阶技术实现了基于低阶不完备模态的模型修正。而对于第2)类问题，尚未展开细致研究。线性“病态”方程组的求解问题在数学领域并不罕见，常用的手段是正则化方法。最早的正则化方法由Tikhonov[7]1963年提出，通过正则化参数的选择，使得求解在最小二乘解和最小范数解之间取得平衡，从而控制解的振荡和方程的拟合程度。目前最常用的正则化方法主要依据矩阵的奇异值分解(SVD)技术，如截断奇异值分解(TSVD)方法[8-10]。文献[11]对各种正则化方法进行了较系统深入的研究，针对中小规模、大规模等不同的“病态”系统提出了适用的正则化方法。研究还指出，利用正则化方法求解不适定问题的关键是正则化参数的选择，这一过程通常需要借助一些数学统计手段，如Morozov偏差准则、广义交叉检验法、“L”曲线法等，这些方法均需要进行大量的试算。正则化方法尽管在理论上比较完善，但在模型修正的实际应用中却仍然存在一些问题，如对于不同问题，同一算法表现出的收敛性却不同，很难用一个通用的方法解决不适定问题[12]。

将这种新的模型修正方法(CMCM 方法)应用到海上风电支撑结构的损伤识别中。首先，给出了该方法的推导过程；然后介绍了正则化求解方法和Guyan扩阶技术，分别用以解决该方法在测量噪声影响下存在的“病态”求解问题和实测空间不完备问题；最后，以一个三脚架式海上风电支撑结构为模型进行数值模拟验证。为提高损伤识别精度，文章结合算例提出一种两步式模型修正过程，通过第一步模型修正实现对损伤的初步定位；在此基础上，进行第二步模型修正，从而实现对损伤位置和损伤程度的精确识别。结果表明，在实测模态不完备且含有测量噪声条件下，结合正则化求解方法和Guyan扩阶技术，采用两步CMCM方法可以准确地识别结构的损伤。

1 基本理论

1.1CMCM方法

将未损伤结构的整体刚度矩阵和质量矩阵记为K和M，一般可通过有限元模型分析得到，其对应于第i阶特征值λi和特征向量Φi的关系可表示为

假设损伤结构模型的刚度矩阵记为K*，并将其写为对未损伤模型刚度矩阵K的一个修正：

其中，Kn为第n个单元在整体坐标系下的刚度矩阵；Ne为待识别损伤单元的个数；αn为刚度修正系数，即损伤程度，可模拟为对应单元的材料弹性模量的变化。另外，假定结构在损伤过程中，不发生质量的变化，即

对于损伤模型，同样存在特征关系：

将式(4)前乘Φi的转置(Φi)T，可得

将式(2)、式(3)代入式(5)得到一线性方程：

用m代替ij，式(6)改写为

由以上推导过程可以得出，当有Ni阶模态从未损伤模型中获取，并且Nj阶模态从损伤结构的模态测试中获得时，共可以组建Nm=Ni×Nj个形如式(7)的线性方程，将其写为方程组

其中，A为Nm×Ne的矩阵；x为包含Ne个刚度损伤系数αn的列向量；b为Nm×1的列向量。

求解模型修正方程组(8)即可实现对结构进行损伤定位和损伤程度的识别。由于未损伤模型的模态是通过有限元分析获得，其可用模态数较多，因此该方法所构建的方程组通常属于超定情况，即Nm≥Ne，此时可采用奇异值分解法来求取方程组的最小二乘解。

对式(8)的系数矩阵A做奇异值分解：

其中，U∈Nm×Nm，V∈Ne×Ne为酉矩阵，Σ∈Nm×Ne为奇异值矩阵。由于矩阵A的某些奇异值可能为0，因此奇异值矩阵Σ可写为

其中，对角阵Σr=diag(σ1,σ2,…,σr)的对角元素为矩阵A的奇异值，且满足σ1≥σ2≥…≥σr>0，r为A的秩。

将酉矩阵U和V分别记为U=(u1,u2,…,uNm)，V=(v1,v2,…,vNe)，其中u1,u2,…,uNm和v1,v2,…,vNe分别为酉矩阵U和V的列向量，则式(10)可写为

方程(8)的解可写为：

1.2正则化求解方法

Tikhonov[7]将最小二乘解和最小范数解综合考虑：

其中，Γ为正则化矩阵，很多情况下取Γ=I或Γ=αI，其中α称为正则化参数，用于控制最小二乘解和最小范数解之间的平衡，当α值较大时，可以保证‖x‖值较小，但可能会使得解不满足方程组Ax=b；反之，如果α值较小，式(13)较接近原问题，虽然能满足方程组的解，但该解会发生振荡，使得没有物理意义。α的取值对求解非常关键。

假定A的奇异值分解为式(9)，则方程组(8)的正则解可写为

称为正则化因子，其作用是过滤掉小奇异值对解的贡献，从而起到稳定解的作用。为方便计算，在截断奇异值分解法(TSVD)中，直接舍弃式(12)中的小奇异值。即取过滤因子为：

如果取σ1≥…≥σk≥α≥σk+1≥…≥σr≥0，则正则解为

其中，k称为截断数，当k=rank(A)时，式(17)即为式(12)所示最小二乘解。

利用正则化方法求解病态方程组的关键是选择正则化参数。对于正则化参数α∈(0,∞)或k∈(1,r)，曲线(ln(‖Γx‖2)，ln(‖b-Ax‖2))形状呈“L”型，称为“L-曲线”，其角点正好使得最小二乘解和最小范数解得到平衡。因此可以通过该曲线的最大曲率来确定角点，从而确定截断数k。

1.3Guyan扩阶

受海上工作环境的限制，在实际测试中，传感器布置数目非常有限，这导致实测模态信息往往是空间不完备的。采用Guyan方法对实测模态进行扩阶处理。

将结构自由振动方程写为

其中，下标m表示主自由度(测试自由度)，s表示从自由度(未测试自由度)。忽略惯性力，有

1.4模态应变能法

为验证方法的优越性，拟采用模态应变能法[13]进行对比研究，基本公式如下：

对于含有ne个单元的线性无阻尼系统，结构的第i阶模态应变能定义为：

第j个单元的第i阶模态应变能分量写为：

当有nm阶实测模态时，第j个单元的损伤指标βj可以写为：

当αj<0时，表示第j个单元有损伤。

2 两步式模型修正过程

当实测模态空间不完备且包含测量噪声时，振型扩阶和正则化求解过程难免会造成求解误差，该误差将随着实测模态的不完备程度和噪声水平的提高而增大，因此，单纯地采用CMCM方法进行一次模型修正难免会造成损伤的误判。通过大量试算，结合CMCM方法，提出一种两步式模型修正过程来精确识别结构损伤。

第一步：将结构的全部单元选为损伤识别(修正)对象，利用CMCM方法构建方程组，根据求解结果，结合物理意义，初步定位损伤位置；

第二步：将初步定位损伤的单元选为损伤识别(修正)对象，利用CMCM方法再次进行识别。根据本次结果通常可以较准确地识别损伤位置和损伤程度。

3 数值算例

3.1算例介绍

一个3.0 MW的海上风电装置，采用三脚架式(tripod)基础结构，如图1所示。应用Matlab编程，建立该风电结构的有限元简化模型，如图2所示，共22个单元，17个结点，杆件采用两结点三维梁单元；结点13和17处分别施加6 t和111 t的集中质量，分别模拟工作平台处设备和塔顶风机的质量。桩基采用等效固定桩来模拟，等效桩底距泥面的长度取为10 m[14-16]。结构材料为钢材，弹性模量取为E=2.1×1011N/m2，密度ρ=7 850 kg/m3，采用集中质量矩阵。

通过特征值分解，得到有限元基准模型的模态参数，如图3所示，其前六阶频率分别为：0.568 5 Hz、0.569 4 Hz、2.787 3 Hz、2.790 2 Hz、5.023 7 Hz和5.034 Hz。

图1 三脚架式海上风电支撑结构Fig.1 Tripod offshore wind turbine structure

图2 结构的有限元简化模型Fig.2 FE model of the structure

图3 有限元模型的前六阶模Fig.3 The first six modes of FE model

3.2损伤模拟

在本研究中，损伤程度定义为某个构件或单元的刚度损失百分比。比如单元13损伤了20%，即该单元的损伤后刚度比基准模型减少了20%。在该算例中，假定3个单元发生了损伤，分别为三脚架一斜撑(单元13)，水面附近立柱(单元17)，塔架靠近工作平台段(单元19)，损伤程度均设为30%。损伤结构是通过修改有限元模型中这3个单元的刚度矩阵来进行模拟的，通过特征值分解，得到损伤模型的模态频率和振型来模拟实测模态参数。损伤模型的前六阶频率分别为：0.520 5 Hz、0.525 5 Hz、2.691 9 Hz、2.733 5 Hz、4.927 9 Hz和4.938 Hz，与损伤前结构相比，频率均有所降低。

3.3损伤识别

3.3.1 工况1：振型空间完备且不含噪声

为验证CMCM方法的适用性，首先假定通过振动测试可以得到损伤结构的前6阶模态，振型空间完备，且不含测量噪声。

将全部22个单元作为损伤识别对象，即Ne=22，利用前6阶实测模态和前6阶有限元模态，可构建Nm=36个CMCM方程。此时方程个数大于未知数个数，为超定方程，采用奇异值分解法(SVD)求取其最小二乘解，结果如图4所示，损伤位置和损伤程度均可实现准确的识别。

3.3.2 工况2：振型空间完备且包含噪声

在实际模态测试中，模态频率的识别结果一般较准确，而振型的识别结果常常受测量噪声的影响较严重，因此在本算例中假定实测模态振型含有噪声。噪声添加方式为Φi,j=Φi,j(1+εRi,j)。其中，ε为噪声水平，Ri,j为均值为0、方差为1的高斯随机数，Φi,j为对应于第i个自由度的第j阶振型值。即在通过数值模拟获得的实测模态振型的基础上添加一定水平的高斯白噪声。

首先假定噪声水平为ε=0.5%。同样将全部22个单元作为损伤识别对象，利用前6阶实测模态和前6阶有限元模态，构建36个CMCM方程。由于该方程组的条件数较大(cond(A)=3.582 1×105)，噪声的存在，使得最SVD求解结果发生振荡，损伤识别结果失效。因此尝试采用TSVD正则化方法进行求解，并用“L-曲线”法选取截断数k。通过大量试算，画出该系统的“L-曲线”，如图5所示。将曲线中位于角点左侧最靠近的点选为截断数，此时k=13。根据这一截断数，利用TSVD方法进行求解，损伤识别结果如图6(a)所示。单元13、单元17和单元19很明显可以诊断到损伤，但与之相连的单元18和单元20也识别到较大的刚度损伤。因此需要结合物理意义定位损伤位置。通常当某个单元含有刚度损伤时，与之相连的单元也会有模态能量的变化，即也会诊断出一定的刚度损伤[16]。根据这一原则，可以排除单元20(与该单元相连的单元21并未诊断出损伤)；但单元18是否含有损伤尚不能确定。因此，初步定位损伤单元为单元13、单元17、单元18和单元19。

图4 工况1损伤识别结果Fig.4 Damage identification result of Case 1

图5 工况2“L-曲线”图Fig.5 L-curve of Case 2

将初步定位的四个损伤单元作为识别对象，利用前6阶实测模态和前6阶有限元模态，重新构建Nm=36个CMCM方程，此时Ne=4。由于此时未知数个数很少，直接采用奇异值分解法(SVD)求取其最小二乘解，结果如图6(b)所示，在0.5%噪声影响下，经过第二步识别，可以相当准确的定位出损伤，且对损伤程度的求解误差也非常小。

图6 利用两步式CMCM模型修正方法识别工况2损伤Fig.6 Damage identification with two-step CMCM model updating method of Case 2

为探讨该两步式方法在不同测量噪声水平下的适用性，将噪声水平分别设为1%、 5%、10%和20%进行对比研究，如图7所示。损伤识别方程组的构建方式和噪声添加方式同上。采用两步式模型修正过程进行损伤识别，结果如图8～10所示。

图7 L-曲线图Fig.7 L-curves

图8 利用两步式CMCM模型修正方法进行损伤识别 (1%噪声)Fig.8 Damage identification with two-step CMCM model updating method (1% noise)

图9 利用两步式CMCM模型修正方法进行损伤识别 (5%噪声)Fig.9 Damage identification with two-step CMCM model updating method (5% noise)

图10 利用两步式CMCM模型修正方法进行损伤识别 (10%噪声)Fig.10 Damage identification with two-step CMCM model updating method (10% noise)

通过比较可以看出，随着噪声水平的提高，“L-曲线”逐渐趋于平滑，当噪声水平达到20%时，很难选定截断数k。但对于能够通过“L-曲线”选定出截断数的三个噪声工况，均可以在第一步模型修正中通过TSVD求解方法定位损伤位置，在第二步模型修正中确定损伤程度。但随着噪声水平的提高，损伤程度识别误差不断增大，造成这一结果的主要原因是，随着噪声水平的提高，模型修正方程组的病态性不断增大，为防止解的“震荡”，利用L曲线选取的截断数不断变小，这使得解趋向于最小范数解，而偏离最小二乘解，从而造成误差的增大。

为体现以上两步式模型修正方法的优越性，采用模态应变能法分别对四种噪声水平工况进行了对比研究，结果如图11所示。可以看出，当噪声水平为0.5%时，模态应变能指标可以准确的定位单元13、17和19的损伤，损伤程度均在-0.3附近，与预设值吻合较好；当噪声水平提高到1%时，损伤程度识别误差增大，且明显高于图8所示两步式模型修正方法给出的误差值；当噪声水平增大到5%和10%时，仅能定位到单元17的损伤，且损伤程度识别误差较大。由此可见，两步式模型修正方法比模态应变能法具有更强的噪声鲁棒性。除此之外，模态应变能指标还造成了单元4、5和6的损伤误判，而两步式模型修正方法通过第一步初步定位可以有效地避免这一问题。

图11 利用模态应变能法识别各噪声水平工况下的损伤Fig.11 Damage identification of different noise levels with modal strain energy method

3.3.3 工况3：振型空间不完备且包含噪声

假定实际测试中，只在结点11～17的三个平动自由度(x、y和z)上布置传感器，即通过测试仅获得21个自由度的振型值。且假定实测模态振型中包含1%的噪声。采用2.3中Guyan方法对实测振型进行扩阶。

第一步：先将全部22个单元作为损伤识别对象，利用前6阶实测模态和前6阶有限元模态，构建36个CMCM方程。利用“L-曲线”选取截断数k=8，如图12(a)所示。采用TSVD 方法求取正则解，根据这一结果可初步定位损伤位置为单元13、单元17、单元18和单元19，如图12(b)所示。

第二步：将初步定位的四个损伤单元作为识别对象，利用第3～6阶实测模态和第3～6阶有限元模态，重新构建Nm=16个CMCM方程，此时Ne=4。采用奇异值分解法(SVD)进行求解，结果如图12(c)所示，经过第二次识别，可以相对准确地定位出损伤位置和损伤程度。与图8所示振型空间完备工况相比，损伤识别误差稍大。由此可见，当实测模态振型空间不完备时，扩阶过程也会造成一定的识别误差，但影响不大。

图12 利用两步式CMCM模型修正方法识别工况3损伤Fig.12 Damage identification of Case 3 with two-step CMCM model updating method

图13 利用模态应变能法识别工况3损伤Fig.13 Damage identification of Case 3 with modal strain energy method

同样采用模态应变能法进行对比研究，结果如图13所示。该方法仅能定位到单元17和单元19的损伤，且损伤程度识别结果误差较大；另外，还造成了水平撑单元7～9和单元16等更多单元的损伤误判。与图11(b)所示振型空间完备工况对比可见，模态扩阶过程对模态应变能方法的损伤识别结果影响较大。

4 结语

1)将一种新提出的交叉模型交叉模态模型修正方法(CMCM)应用到了海上风电支撑结构的损伤识别中。并提出了一种两步式模型修正过程：第一步，将结构的所有单元作为对象，利用CMCM方法进行伤识别，根据结果对损伤位置进行初步定位；第二步，将初步定位损伤的单元重新选作损伤识别对象，再次利用CMCM方法实现对损伤位置和损伤程度的精确识别。数值算例表明，该方法与模态应变能法相比具有较强的噪声鲁棒性和振型扩阶容差性，且能有效的避免损伤的误判。

2)当用CMCM方法构建的损伤识别方程组呈现一定程度的“病态”时，在测量噪声的影响下，最小二乘解容易导致求解失效。采用TSVD正则化方法进行求解可有效解决这一问题。

3)在测量模态不完备且含有一定噪声的情况下，结合传统的Guyan扩阶技术，采用两步式模型修正过程可以较好地识别海上风电支撑结构的损伤。

[1] LIU F,LI H,LI W,et al.Experimental study of improved modal strain energy method for damage localization in jacket-type offshore wind turbines [J].Renewable Energy,2014,72:174-181.

[2] SOHN H,FARRAR C R,HEMEZ F M.A review of structural health monitoring literature:1996-2001[R].Los Alamos National Laboratory Rep.No.LA-13976-MS,2003.

[3] PEETERS B,MAECK J,ROECK G D.Vibration-based damage detection in civil engineering:excitation sources and temperature effects [J].Journal of Smart Materials and Structures,2001,10(3):518-527.

[4] HU S J,LI H J,Wang S Q.Cross-model cross-mode method for model updating [J].Mechanical Systems and Signal Processing,2007,21(4):1690-1703.

[5] LI H J,WANG J R,HU S J.Using incomplete modal data for damage detection in offshore jacket structures [J].Ocean Engineering,2008,35(17-18):1793-1799.

[6] 李英超.基于模态参数识别的海洋平台结构模型修正技术研究[D].青岛：中国海洋大学，2012.(LI Yingchao.Study on finite element model updating for offshore platform structures based on vibration test [D].Qingdao:Ocean University of China,2012.(in Chinese))

[7] TIKHONOV A N.Solution of incorrectly formulated problems and the regularization method [J].Soviet Math,Dokl,1963,4：1035-1038.

[8] MARES C,FRISWELL M I.MOTTERSHEAD J E.Model updating using robust estimation [J].Mechanical Systems and Signal Processing,2002,16(1)：169-183.

[9] FRISWELL M I,MOTTERSHEAD J,AHMADIAN H.Finite-element model updating using experimental test data:parameterization and regu1arization[J].Philosophical Transactions of the Royal Society of London,2001,359：169-186.

[10] CALVETTI D,REICHEL L,SGALLARI F,et al.A regularizing lanczos iteration method for underdetermined linear systems [J].Journal of Computational and Applied Mathematics,2000,115：101-120.

[11] 吴颉尔.正则化方法及其在模型修正中的应用[D].南京：南京航空航天大学,2007.(WU Jieer.The regularization method with application to finite element model updating [D].Nanjing:Nanjing University of Aeronautics and Astronautics,2007.(in Chinese))

[12] 肖庭延,于慎根,王彦飞.反问题的数值解法[M].北京：科学出版社,2003.( XIAO Tingyan,YU Shengen,WANG Yanfei.Numerical solution of inverse problem [M].Beijing:Science Press,2003.(in Chinese))

[13] KIM J T,STUBBS N.Improved damage identification method based on modal information[J].Journal of Sound and Vibration,2002,252(2):223-238.

[14] ZAAIJER M B.Foundation modeling to assess dynamic behaviour of 0ffshore wind turbines[J].Applied Ocean Research,2006,28,45-57.

[15] JAN V T.Design of supports for offshore wind turbines [D].Delft:Delft University of Technology,2006.

[16] WANG S Q.Damage detection in offshore platform structures from limited modal data [J].Applied Ocean Research,2013,41:48-56.

Damage identification of offshore wind turbine support structure with a two-step model updating process

LI Yingchao1,ZHANG Min2,WANG Shuqing2

(1.College of Civil Engineering,Ludong University,Yantai 264025,China; 2.College of Engineering,Ocean University of China,Qingdao 266100,China)

P752; P741

10.16483/j.issn.1005-9865.2016.03.004

1005-9865(2016)03-0028-10

20105-06-09

国家自然科学基金项目(51379196，51209189)；山东省自然科学基金项目(ZR2013EEQ006)；泰山学者工程专项经费；鲁东大学科研基金项目(LY2013027)