基于子结构技术和Landweber迭代正则化的载荷识别

2021-08-13周凤缪炳荣李泽杨昌休

机械制造与自动化 2021年4期

周凤，缪炳荣，李泽，杨昌休

(1. 重庆交通职业学院，重庆402247；2. 西南交通大学牵引动力国家重点实验室，四川成都 610031)

0 引言

载荷识别在机械动态结构的健康监控、可靠性分析和故障诊断中都起着至关重要的作用[1]。然而，在许多工程实践中，例如风机叶片上的风力载荷、飞机上的空气动载荷、铁路轨道和车轮之间的相互作用力等，这些载荷由于技术或经济的限制直接进行测量是很困难的，甚至是不可能的。所以，利用载荷识别技术来反求载荷的时间历程和位置是非常有必要的。

在过去的几十年里，动态载荷识别问题已经得到了广泛的研究。提出的方法中主要包含频域方法[2-7]和时域方法[8-13]两种。与频域方法相比，时域方法可以识别包括冲击载荷在内的各种载荷，识别结果具有明确的物理意义。时域方法在工程中具有良好的应用前景。

载荷识别作为结构动力学中第二类逆问题，具有典型的不适定性，无论是对于频域方法还是时域方法，都需要正则化技术来找到准确的解决方案。目前，常用的正则化方法有截断奇异值分解(TSVD)、Tikhonov正则化和迭代正则化方法[14-17]。

上述的载荷识别方法基本上都是关于简单的梁和板上所作用载荷的识别。但关于复杂大型结构的载荷识别研究较少。本文在文献[18]的基础上，利用Landweber迭代正则化求解反问题，以提高载荷识别的精度以及抗噪性，并利用模型减缩技术中的子结构分析，对大型有限元结构进行降阶处理，得到其子结构模型，在保证响应计算精度的同时，极大地缩短求解时间，提高计算效率。

另外，目前载荷识别在铁道车辆中的运用主要集中在轮轨力上的识别[19]，而关于车体表面动态载荷的识别研究还基本处于空白。本文利用提出的子结构技术与Landweber迭代正则化来识别比例车体有限元模型上的载荷。

1 载荷识别原理

1.1 正问题模型的建立

具有n自由度的比例阻尼系统的动力学方程为

(1)

假定系统对单位脉冲的响应，即由载荷作用点到响应测量点的Green函数G(t)，则根据叠加原理，系统所受外载荷所引起的动态响应为

(2)

考虑零初始条件系统，将上式卷积积分以Δ(t)为离散的采样时间间隔，m为采样点数，离散为一组线性方程组，其矩阵形式表示如下：

(3)

进一步，式(3)可以简单表示为

Y=GF

(4)

其中：yi、Gi、fi分别为t=i△(t)的响应、Green函数和待求的外载荷；G为下三角矩阵，是由有限元求解得到的Green脉冲核函数响应建立的。

当响应中含有噪声时，通过方程式(4)求解外部载荷就变成一个病态的、不适定问题，就需要正则化技术解决这类病态问题。

1.2 Landweber迭代正则化方法

考虑抽象的算子方程：

Ax=y

(5)

其中A:X→Y为有界线性算子，X和Y均为可分的Hilbert空间。

Landweber迭代法的迭代格式为

xk=xk-1+wA*(y-Axk-1)

(6)

fk=fk-1+wGT(y-Gfk-1)

(7)

式中的k与式(6)中的一样为迭代次数，发挥着正则化参数的作用。初始解向量f0可设置为f0=0。

1.3 子结构分析方法

在利用上述迭代方法进行载荷识别时，需要计算有限元仿真系统的动响应。对于大型复杂的结构，其有限元模型的网格质量决定了仿真的精度，而越密的网格会导致求解动响应的时间也就越长，有时候甚至求解困难。采用子结构法将大型复杂结构按照一定的原则划分成若干个子结构，先将每一个子结构减缩为很少自由度的模型，然后根据各子结构之间关系，组装得到整体结构的动态特性。本文采用的子结构分析方法属于模型减缩技术中的一种。

子结构法就是将一组单元用矩阵凝聚为一个超单元的过程。系统的运动平衡方程为

(8)

将u分为主自由度us和从自由度uc两部分，并令

u=Tus

(9)

其中T为转换矩阵。

采用 Guyan缩减方法，将相同的转换矩阵T扩展到式(8)中的质量矩阵和阻尼矩阵，即

Ms=TTMT,Cs=TTCT

(10)

新的平衡方程变为

(11)

此时自由度数缩减为主自由度数。

2 数值算例

为了验证上述Landweber迭代正则化识别的准确性、抗干扰性以及子结构分析技术对识别过程中计算效率改善的正确性，选择大型复杂的比例车体结构作为数值仿真模型。在ANSYS中建立如图1所示的有限元模型，主体结构采用shell63板壳单元进行离散化网格划分，该比例车体模型是相对高速列车1∶8的简化模型，长×宽×高为3 022 mm×407 mm×350.8 mm，模型材料为铝合金，即弹性模量E=6.9×1010Pa，密度ρ=2 700 kg/m3，泊松比v=0.3。

图1 比例车体有限元模型

2.1 子结构的应用及结果分析

ANSYS软件中应用Guyan缩减法来计算缩减矩阵，而选取主自由度在矩阵缩减分析中起到关键作用，其准确程度将取决于主自由度的数量和分布位置。因此，根据文献[20]中主自由度的选取方法对该比例车体结构有限元模型进行主自由度的选取。在模型上选取包括加载点和响应测量点在内的共301个节点，缩减自由度后的比例车体有限元模型如图2所示。分析缩减前后结构的模态频率以及响应计算效率。

图2 缩减自由度后的比例车体模型

1)缩减前后模态对比

表1为比例车体缩减前全自由度和缩减后模态频率的对比情况(去掉前6阶刚体模态)。图3为缩减前后模态频率对比图。从表1和图3都可以看出，所选取主自由度后的子结构模型能满足模型的准确度，且对于前4阶模态，选取前后模态频率的相对偏差基本低于1%或在1%附近。这是由于主自由度的选取数量对低价模态影响较小，而对高阶模态影响较大。

表1 比例车体主自由度选取前后模态对比

图3 比例车体模态频率对比图

2)核函数响应计算效率分析

分别对比例车体原模型和子结构缩减模型求脉冲核函数响应，比较这两种模型中计算相同时间历程的核函数响应所需要的实际CPU时间，如表2所示。其中采样周期为Δt=0.000 5s。

表2 缩减前后模型核函数响应计算时间

单位：s

由表2可见，相比于全自由度的原模型，子结构模型大大减少了核函数响应计算时间，且所需计算的响应时间历程越长，子结构缩减模型提高计算效率的效果就越明显。因此，应用子结构技术将大大提高载荷识别过程中计算核函数响应的效率。

2.2 载荷识别结果分析

1)不同形式载荷下的识别结果

在比例车体子结构模型的节点5614上施加垂直于上顶板的载荷，即z向上加载，并选取距离测量点较近的响应节点5795(z向)上的位移作为响应用于载荷识别。所施加的载荷分别取如下几种形式：

1)正弦载荷：

(12)

2)三角载荷：

(13)

3)方波载荷：

(14)

正弦载荷、三角载荷和方波载荷选用不同的采样频率，即f1=2000Hz、f2=4000Hz和f3=3 333.33Hz，利用Landweber迭代正则化方法进行载荷识别，并对比传统的Tikhonov正则化。各载荷形式下的识别结果分别如图4-图6所示。

由图4-图6可以看出，对于不同形式的载荷，在没有噪声的情况下，传统的Tikhonov正则化方法和Landweber迭代正则化都能准确地识别载荷，这是由于在响应数据不含有噪声时(实际情况中几乎是不存在的)，可以直接通过矩阵求逆法反求载荷，即此时有无正则化的识别效果相差不大，所以这两种正则化方法下的载荷识别结果相似。但从图4和图6也可以发现，传统Tikhonov正则化方法在载荷值突变处以及时间起始点处，相比于Landweber迭代正则化的识别误差较大。

图4 正弦载荷识别结果对比

图5 三角载荷识别结果对比

图6 方波载荷识别结果对比

2)不同噪声水平下的识别结果

由于实际测量的响应信号y中总是包含噪声，在仿真得到的响应数据中加入不同水平的随机噪声来模拟测试误差。带噪声的响应可表示成如下形式：

ynoise=y(t)+lnoisestd[y(t)]rand(-1,1)

(15)

其中：y(t)为仿真得到的位移响应；std[y(t)]为位移响应y(t)的标准差；lnoise为噪声水平的百分数；rand(-1,1)为区间[-1,1]上的随机数。

首先，对于正弦载荷在噪声水平分别为3%、5%和10%的载荷识别结果如图7所示。从图7可以看出，随着噪声水平的增大，Tikhonov正则化识别的误差逐渐增加，并且到噪声水平为10%时，已经不能识别出所加载荷。图7(c)中，Tikhonov正则化方法在t=0.03s附近，识别结果已经趋于发散。而Landweber迭代正则化在这3种噪声水平下都能较准确地识别，识别精度并没有随着噪声水平的波动而发生较大变化。

图7 在不同噪声水平(3%、5%、10%)下正弦载荷的识别结果

为了更进一步说明噪声的影响，选取正弦载荷在噪声水平为5%下，随机5个不同的离散时间点上识别的结果如表3所示，得到Landweber迭代正则化识别误差，除了少许的几个时间点，基本在5%以下，且大部分的识别相对误差低于Tikhonov方法。同时分析在噪声水平为5%下，三角载荷和方波载荷的识别结果如图8、图9所示。从这两幅图中可以看出，与正弦载荷类似，在5%的噪声水平下Landweber正则化识别三角和方波这两种载荷的识别精度明显高于Tikhonov正则化。