CN117669211A - 基于参数化时域传递率的结构参数辨识及方差计算方法 - Google Patents
基于参数化时域传递率的结构参数辨识及方差计算方法 Download PDFInfo
- Publication number
- CN117669211A CN117669211A CN202311664658.8A CN202311664658A CN117669211A CN 117669211 A CN117669211 A CN 117669211A CN 202311664658 A CN202311664658 A CN 202311664658A CN 117669211 A CN117669211 A CN 117669211A
- Authority
- CN
- China
- Prior art keywords
- matrix
- formula
- calculating
- variance
- index
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Pending
Links
- 238000004364 calculation method Methods 0.000 title claims abstract description 21
- 239000011159 matrix material Substances 0.000 claims abstract description 173
- 238000000034 method Methods 0.000 claims abstract description 62
- 230000004044 response Effects 0.000 claims abstract description 59
- 238000013016 damping Methods 0.000 claims abstract description 42
- 230000036541 health Effects 0.000 claims abstract description 12
- 238000012544 monitoring process Methods 0.000 claims abstract description 12
- 239000013598 vector Substances 0.000 claims description 24
- 238000005314 correlation function Methods 0.000 claims description 19
- 230000001133 acceleration Effects 0.000 claims description 12
- 125000003275 alpha amino acid group Chemical group 0.000 claims description 12
- 150000001875 compounds Chemical class 0.000 claims description 12
- 230000035945 sensitivity Effects 0.000 claims description 10
- 230000000903 blocking effect Effects 0.000 claims description 8
- 238000000638 solvent extraction Methods 0.000 claims description 6
- 238000006243 chemical reaction Methods 0.000 claims description 4
- 238000006073 displacement reaction Methods 0.000 claims description 4
- 238000005070 sampling Methods 0.000 claims description 4
- 238000010276 construction Methods 0.000 claims description 3
- NAWXUBYGYWOOIX-SFHVURJKSA-N (2s)-2-[[4-[2-(2,4-diaminoquinazolin-6-yl)ethyl]benzoyl]amino]-4-methylidenepentanedioic acid Chemical compound C1=CC2=NC(N)=NC(N)=C2C=C1CCC1=CC=C(C(=O)N[C@@H](CC(=C)C(O)=O)C(O)=O)C=C1 NAWXUBYGYWOOIX-SFHVURJKSA-N 0.000 claims description 2
- 230000017105 transposition Effects 0.000 claims description 2
- 238000005516 engineering process Methods 0.000 abstract description 16
- 238000005259 measurement Methods 0.000 abstract description 9
- 230000003595 spectral effect Effects 0.000 description 9
- 238000011156 evaluation Methods 0.000 description 7
- 238000004458 analytical method Methods 0.000 description 5
- 238000000354 decomposition reaction Methods 0.000 description 4
- 238000012937 correction Methods 0.000 description 2
- 230000000694 effects Effects 0.000 description 2
- 238000012423 maintenance Methods 0.000 description 2
- 239000000463 material Substances 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 238000012546 transfer Methods 0.000 description 2
- 230000015556 catabolic process Effects 0.000 description 1
- 238000006731 degradation reaction Methods 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 238000010606 normalization Methods 0.000 description 1
- 230000003449 preventive effect Effects 0.000 description 1
- 238000011160 research Methods 0.000 description 1
Classifications
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
- Y02T90/00—Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation
Landscapes
- Measurement Of Mechanical Vibrations Or Ultrasonic Waves (AREA)
Abstract
本发明公开了基于参数化时域传递率的结构参数辨识及方差计算方法,具体涉及结构健康监测技术领域,首先测量两个工况下结构振动响应并划分为参考信号和非参考信号,在两个工况下构建参数化模型,并计算参数化模型的系数矩阵及其方差;根据系数矩阵构建友矩阵,对友矩阵进行特征值分解并据此计算结构模态频率、阻尼比和模态振型;通过系数矩阵的方差,计算友矩阵方差。本发明与已有的直接利用振动响应信号进行结构参数辨识的技术相比,受到结构振动响应测量噪声的影响更小,结构参数辨识精度和可靠性更高。
Description
技术领域
本发明涉及结构健康监测技术领域,具体涉及基于参数化时域传递率的结构参数辨识及方差计算方法。
背景技术
我国在交通、楼宇、航空航天等行业的基础设施建设井喷式发展,海量工程结构的安全保障和维护管理成为新的难题。通过测量结构在运行状态下的振动响应数据,从振动响应数据中辨识结构参数,从而开展工程结构健康监测、状态评估和智能控制,是尽早发现工程结构损伤、延缓结构性能退化的重要方法,对于实现工程结构预防性维护管理、减少财产损失具有重要意义。
从测量的结构振动响应数据中辨识结构参数的技术称为工作模态分析技术,主要辨识的结构参数包括模态频率、阻尼比和模态振型。该技术目前已成为结构健康监测领域的研究热点之一,涌现了一大批结构参数辨识算法。工作模态分析技术仅利用结构在运行状态下的振动响应数据,而假设作用在结构的载荷为白噪声序列,因此当结构的载荷不为白噪声序列时,容易导致虚假结果,给结构健康监测和状态评估带来误判,严重时会引起重大的经济损失。
目前,部分已有技术可解决非白噪声序列载荷下的结构参数辨识问题。其中,基于响应传递率和功率谱密度传递率的方法在工程中应用广泛。响应传递率在结构模态频率处只与结构模态振型有关,通过联合多个工况下的响应传递率函数可辨识结构参数。在已有技术中,基于响应传递率的方法可分为时域响应传递率方法和频域传递率方法,其中时域响应传递率方法的精度更高。基于响应功率谱密度的方法通过联合多个转换通道下的响应功率谱密度函数,辨识结构参数。
然而,实际工程应用表明,基于响应传递率的方法容易受到振动响应数据中测量噪声的影响,而基于功率谱密度传递率的方法在很多情况下容易失效。此外,目前基于响应传递率的方法和基于功率谱密度传递率的方法都无法计算结构参数的不确定度。
发明公开号CN107561934 A公开了一种基于多参考时域传递率的仅输出结构模态参数辨识方法,该方法将采集的工程结构响应信号分成非参考输出信号与参考输出信号两组,分别作为多参考时域传递率模型的输入和输出,能处理非白输入条件下的结构模态辨识问题,结构不受载荷形式限制。该方法虽然相对于传统的频域方法具有优势,但仍然存在辨识结果受振动响应测量噪声影响大、无法计算结构参数方差的问题,而结构参数的不确定度是评估结构参数辨识效果和衡量结构参数受噪声干扰程度的重要特征,并且结构参数的不确定度在结构模型修正、健康监测及状态评估中具有重要的应用。在实际工程中,结构参数的不确定度可通过结构参数的方差进行表征,方差越小,表明结构参数辨识效果越好,辨识结果可靠性越高,结构参数受噪声干扰程度越低。
本发明所公开的技术是为了解决发明公开号CN107561934 A中存在的辨识结果受振动响应测量噪声影响大、无法计算结构参数方差的问题而提出的。
发明内容
为此,本发明提供基于参数化时域传递率的结构参数辨识及方差计算方法,以解决背景技术中提出的问题。
本发明的目的是通过下述技术方案实现的:
本发明公开的基于参数化传递率的结构参数辨识及方差计算方法,首先测量两个工况下结构振动响应并计算相关函数,随后将相关函数划分为参考信号和非参考信号,在两个工况下构建参数化模型,并计算参数化模型的系数矩阵及其方差;随后,根据系数矩阵构建友矩阵,对友矩阵进行特征值分解并据此计算结构模态频率、阻尼比和模态振型;之后,通过系数矩阵的方差,计算友矩阵方差;接着,计算友矩阵特征值、特征向量相对于友矩阵的灵敏度矩阵,模态频率、阻尼比相对于友矩阵特征值的灵敏度矩阵,以及模态振型相对于友矩阵特征向量的灵敏度矩阵;最后,根据上述计算的灵敏度矩阵及友矩阵方差,计算结构模态频率、阻尼比和模态振型的方差。因此,本发明所公开方法辨识的结构参数包括模态频率、阻尼比和模态振型,并计算出模态频率、阻尼比和模态振型的方差。
具体包括以下步骤:
步骤1:采用位移传感器、速度传感器或加速度传感器,测量两个不同工况下工程结构的位移、速度或加速度振动响应信号;
步骤2:确定转换通道,计算结构振动响应信号的相关函数;
步骤3:确定相关函数中的参考信号和非参考信号,在两个工况下构建参数化模型,估计出参数化模型的系数矩阵;
步骤4:计算系数矩阵的方差矩阵;
步骤5:根据步骤3得到的系数矩阵,构造友矩阵,再计算结构模态频率、阻尼比和模态振型;
步骤6:计算友矩阵的方差矩阵;
步骤7:计算结构模态频率、阻尼比和模态振型的方差。
优选的,步骤1中,采用加速度振动响应信号,两个工况下测量的振动响应信号分别记为x(1)[t]和x(2)[t],其中x(1)[t]和x(2)[t]均是维度为No×1的向量,上标(1)和(2)表示第一个工况和第二个工况,[t]表示第t个时刻,No表示传感器个数并且为大于等于2的整数。
优选的,步骤2的具体包括:
步骤2.1:根据步骤1得到的振动响应信号x(1)[t],选取x(1)[t]中均方根最大的元素,记为其中z为整数且1≤z≤No,表示振动响应信号向量x(1)[t]的第z个元素,则转换通道为z;
步骤2.2:采用公知的相关函数估计方法,分别计算振动响应x(1)[t]与振动响应x(2)[t]与/>的相关函数,分别记为s(1)[t]和s(2)[t],其中,s(1)[t]和s(2)[t]均为维度是No×1的向量。
优选的,步骤3的具体包括:
步骤3.1:在No个传感器中挑选出Nr个传感器,其中Nr为整数并且 表示向下取整;将挑选出的Nr个传感器作为参考传感器,对应的相关函数为参考信号,两个工况下的参考信号分别记为/>和/>其余的Nu=No-Nr个传感器为非参考传感器,对应的相关函数为非参考信号,分别记为/>和/>
步骤3.2:根据步骤3.1得到的和/>构建参数化模型如式(1)所示:
式(1)中,上标“(k)”表示第k个工况,和/>表示第k个工况下的系数矩阵,e(k)[t]表示第k个工况下的残差序列,n表示模型阶数,模型阶数n可利用公知的赤池信息准则(即AIC准则)进行确定;
步骤3.3:式(1)中,所有的系数矩阵和/>由式(2)计算得到
式(2)中,上标“T”和“-1”分别表示矩阵转置和求逆,N为时刻点个数,Θ(k)和Ψ[t]分别为
通过式(2)计算得到Θ(k)后,即可按照式(3)的第一式得到系数矩阵和/> 和/>的矩阵维度分别为Nu×Nu和Nu×Nr。
优选的,步骤4的具体包括:
步骤4.1:根据式(2)计算得到的Θ(k),计算两个工况下的残差序列e(k)[t],如式(4)所示
步骤4.2:根据式(4)计算得到的e(k)[t],计算两个工况下残差序列e(k)[t]的方差矩阵,如式(5)所示
步骤4.3:计算两个工况下系数矩阵Θ(k)的方差矩阵,如式(6)所示
优选的,步骤5的具体包括:
步骤5.1:令为维度Nu×Nu的单位矩阵,即/>其中/>表示维度Nu×Nu的单位矩阵,构造矩阵Ci如式(7)所示
进一步,计算其中上标“+”表示求矩阵伪逆,随后构造友矩阵Γ如式(8)所示
式(8)中,和/>分别表示维度为No×No的单位矩阵和全零矩阵;
步骤5.2:对友矩阵Γ进行特征值分解,记得到的第l个特征值及对应的特征向量分别为λl和则第l阶模态频率fl和阻尼比ξl分别为
式(9)中,Ts为振动响应信号采样的时间间隔,ln(·)和|·|分别表示取自然对数和求绝对值运算;
步骤5.3:根据特征向量计算结构模态振型向量/>首先取特征向量/>的前No个元素,记为/>随后按式(10)计算模态振型向量
式(10)中,为向量/>的第m个元素,Re(·)表示取实部,m的取值是任意的,只需满足/>即可,优选/>值较大的元素;
由步骤5.2和步骤5.3计算得到的模态频率、阻尼比和模态振型即为最终辨识的结构参数。
优选的,步骤6的具体包括:
步骤6.1:将步骤4得到的系数矩阵Θ(k)的方差矩阵进行分块,如式(11)所示
式(11)中,的维度为/> 的维度为(n+1)NuNr×(n+1)NuNr,/>与/>的维度分别为/>和/>根据分块结果,组装得到新的矩阵/>如式(12)所示
步骤6.2:定义两个指标矩阵Eindex和E′index,如式(13)所示
进一步组装得到指标矩阵Eindex,如式(14)所示
步骤6.3:根据步骤6.2得到的指标矩阵Eindex,对其进行分块,如式(15)所示
式(15)中,每个矩阵的维度均为2Nu×Nu,每个/>矩阵的维度均为2Nu×Nr,根据分块结果,对指标矩阵Eindex进行重组,得到新的指标矩阵/>如式(16)所示
步骤6.4:将步骤6.2和6.3得到的指标矩阵Eindex和拉直,分别记为eindex和生成置换矩阵P1,其维度为2(n+1)NuNo×2(n+1)NuNo;生成置换矩阵P1的具体方法为:取出eindex的第ε个元素,记为eε,取出/>的第ε个元素记为/>矩阵P1第/>行第eε列的元素为1,整数ε从1到2(n+1)NuNo依次递增,矩阵P1其余元素均为0;
步骤6.5:计算C0到Cn的方差ΣD,如式(17)所示
对矩阵ΣD进行分块,如式(18)所示
式(18)中,Σ1矩阵的维度为2nNuNo×2nNuNo,Σ4矩阵的维度为2NuNo×2NuNo;
计算矩阵C0′到Cn′的方差ΣD′,如式(19)所示
ΣD′=(MT+Q)Σ4(MT+Q)T+EΣ1ET+(MT+Q)Σ3ET+EΣ2(MT+Q)T (19)
式(19)中,其中D=[C0 C1 … Cn-1];T矩阵满足/>其中vec(·)表示矩阵拉直运算;
步骤6.6:生成维度为(nNo)2×(nNo)2的置换矩阵P2,具体步骤如下:
步骤6.6.1:生成指标矩阵如式(20)所示
式中表示维度为(n-1)No×nNo的零矩阵;将/>拉直,记为/>
取出的第ε个元素,记为/>若/>则置换矩阵P2的第ε行均为0,反之,若则置换矩阵P2第ε行第/>列元素值为1;遍历整数ε从1到2(n+1)NuNo依次递增,矩阵P2其余元素均为0;
计算友矩阵Γ的方差矩阵ΣΓ,如式(21)所示
优选的,步骤7的具体包括:
步骤7.1:计算特征值λl和特征向量相对于友矩阵Γ的灵敏度矩阵/>和/>如式(22)所示
式(22)中,ψl为Γ的第l阶左特征向量,即满足
步骤7.2:计算模态频率fl和阻尼比ξl相对于特征值λl的灵敏度矩阵如式(23)所示
式(23)中,Im(·)表示取虚部;
步骤7.3:计算模态振型对于特征向量/>的灵敏度矩阵,如式(24)所示
式(24)中,通过取式(22)得到的/>矩阵的前No行得到,/>和/>分别表示维度为No′(m-1)和No′(No-m)的零矩阵;
步骤7.4:计算模态频率fl、阻尼比ξl和模态振型的方差,如式(25)所示
式(25)中,取矩阵第一行第一列的元素,即为模态频率fl的方差,取/>第二行第二列的元素,即为阻尼比ξl的方差,/>即为模态振型/>的方差。
优选的,步骤7计算得到的模态频率fl、阻尼比ξl和模态振型的方差可用于评估步骤5中辨识的模态频率fl、阻尼比ξl和模态振型/>的精度。
优选的,辨识的结构模态频率fl、阻尼比ξl和模态振型及其方差应用于工程结构的健康监测及状态评估。
本发明具有如下优点:
1、本发明公开的基于参数化传递率的结构参数辨识及方差计算方法,方法辨识的结构参数不受非白噪声载荷的影响,辨识的结构参数中不包含虚假模态,因此拓展了已有工作模态分析技术的工程应用范围。
2、本发明公开的基于参数化传递率的结构参数辨识及方差计算方法,在辨识结构参数的同时,还能计算得到结构参数的方差,通过结构参数方差对辨识结构参数的精度进行评估。结构参数方差越小,说明辨识精度越高,受到振动响应测量噪声影响越小。因此,本发明计算的结构参数方差提供了重要的辨识精度评估依据,工程实用性更强。
3、本发明公开的基于参数化传递率的结构参数辨识及方差计算方法,在辨识得到结构参数的同时,还能计算得到结构参数的方差,基于辨识的结构参数及其方差,可将其应用于考虑结构参数不确定性的结构健康监测、模型修正及状态评估中,拓展了结构参数辨识技术的工程应用范围,工程可靠性和实用性更强。
4、本发明公开的基于参数化传递率的结构参数辨识及方差计算方法,结合了两个工况下的相关函数对结构参数进行辨识,与已有的直接利用振动响应信号进行结构参数辨识的技术相比,受到结构振动响应测量噪声的影响更小,结构参数辨识精度和可靠性更高。
附图说明
图1为本发明公开的基于参数化传递率的结构参数辨识及方差计算方法的流程图;
图2为具体实施方式中的桁架结构示意图,其中数字表示各杆的编号,结点由带圈数字标出,结点①、②、③、④处的黑点表示集中质量,图右侧的字母至表示各杆的横截面积,单位为cm2;
图3为具体实施方式中两个载荷工况下的桁架结构加速度响应信号的功率谱密度曲线,其中实线为功率谱密度曲线,虚线圆圈表示非白噪声载荷的频率。
图4为具体实施方式中本发明所公开方法辨识出的桁架结构模态振型与置信区间,其中黑色实线为真实的模态振型,黑色点线为辨识出的模态振型,黑色虚线为辨识的模态振型的3σ界。
具体实施方式
以下由特定的具体实施例说明本发明的实施方式,熟悉此技术的人士可由本说明书所揭露的内容轻易地了解本发明的其他优点及功效,显然,所描述的实施例是本发明一部分实施例,而不是全部的实施例。基于本发明中的实施例,本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例,都属于本发明保护的范围。
为了更好地说明本发明的目的和优点,通过对一个非白环境载荷下的桁架结构进行结构模态参数识别,对本发明做出详细解释。
实施例1:
本实施例的桁架结构如图2所示,桁架的尺寸及各杆的横截面积在图中给出。桁架材料的杨氏模量69.8Gpa,材料密度2770kg/m3,在结点1、2、3、4上各附加454kg的集中质量,桁架结构具有八阶模态,八阶结构模态频率分别为5.90、14.98、17.15、20.80、26.81、31.52、43.68、47.54,单位为Hz,阻尼比分别为1.53%、1.00%、1.00%、1.04%、1.14%、1.24%、1.55%、1.66%。桁架结构共含有①、②、③、④共4个未约束结点,每个结点具有沿x轴和y轴的两个自由度,共8个自由度。按结点编号、先x轴后y轴的顺序对桁架结构的自由度进行编号,例如自由度3和4分别表示结点②沿x轴和沿y轴的平移自由度。
考虑两个载荷工况,在第一个载荷工况中,自由度1、2、3和4受到非白噪声载荷的作用,其余4个自由度不受载荷;在第二个载荷工况中,自由度5、6、7、8受到非白噪声载荷的作用,而其余自由度不受载荷作用。在两个载荷工况中,非白噪声载荷均为主频10Hz、阻尼比3%的有色噪声。在两个载荷工况下,桁架结构的加速度响应由Newmark-β方法计算得到,响应信号采样率为256Hz,时间长度为256s,则信号采样时间间隔Ts=1/256s,信号时刻点数目N=65536。在本实施例中,计算得到的加速度响应作为测量的结构振动加速度响应信号,进行桁架结构参数的辨识和方差计算。
两个载荷工况下的桁架结构加速度响应信号的功率谱密度曲线如图3所示。可以看出,两个载荷工况下的响应功率谱密度曲线均包含9个峰值,其中10Hz处的峰值是由非白噪声载荷的主频引起的,而其他8个峰值是由结构模态频率引起的。若采用已有技术中假设结构载荷为白噪声序列的方法对结构参数进行辨识,将会把10Hz辨识成结构模态频率,导致虚假模态,影响结构健康监测与状态评估的可靠性。
本实施例公开的基于参数化传递率的结构参数辨识及方差计算方法,包括以下步骤:
步骤1:在本实施例中,采用桁架结构的加速度响应信号进行结构参数辨识与方差计算。因此,两个工况下振动响应信号x(1)[t]和x(2)[t]均是维度为8×1的向量。
步骤2:确定转换通道z,计算结构振动响应信号的相关函数。
步骤2.1:选取x(1)[t]中均方根最大的元素,在本实施例中为因此转换通道为2。
步骤2.2:分别计算振动响应x(1)[t]与振动响应x(2)[t]与/>的相关函数s(1)[t]和s(2)[t],均为8×1的向量。
步骤3:确定相关函数中的参考信号和非参考信号,在两个工况下构建参数化模型,估计出参数化模型的系数矩阵。
步骤3.1:本实施例中,No=8。挑选自由度1、2、3、4对应的相关函数作为参考信号,其余4个自由度对应的相关函数为非参考信号,因此Nu=Nr=4。参考信号非参考信号/>和/>的维度均为4×1的向量。
步骤3.2:在两个载荷工况下,构建参数化模型如式(1)所示,根据赤池信息准则确定阶数n为15,响应信号的时间长度。
步骤3.3:按照式(3)组装得到矩阵Ψ(1)[t]和Ψ(2)[t],根据式(2)计算得到系数矩阵Θ(1)和Θ(2),均为4×124的矩阵,而系数矩阵和/>的矩阵维度均为4×4。
步骤4:计算系数矩阵Θ(1)和Θ(2)的方差矩阵。
步骤4.1:根据Θ(1)、Θ(2)和式(4)计算两个工况下的残差序列e(1)[t]和e(2)[t]。
步骤4.2:根据e(1)[t]、e(2)[t]和式(5),计算其方差矩阵和/>均为4×4的矩阵。
步骤4.3:根据式(6)计算Θ(1)、Θ(2)的方差矩阵和/>维度均为496×496。
步骤5:计算结构模态频率、阻尼比和模态振型。
步骤5.1:令按式(7)构造矩阵C0到C15。计算/>其中i=0,1,…,14,随后按式(8)构造友矩阵Γ,维度为120×120。
步骤5.2:对友矩阵Γ进行特征值分解,根据得到的特征值和式(9)计算得到模态频率和阻尼比。
步骤5.3:取各特征向量前8个元素作为向量选取第一个元素进行向量归一化,即式(10)中m=1,按式(10)计算模态振型向量。
经过步骤5.2和步骤5.3,辨识得到结构模态频率、阻尼比和模态振型。
步骤6:计算友矩阵Γ的方差矩阵ΣΓ。
步骤6.1:按式(11),将矩阵进行分块,/>矩阵的维度分别为240×240、240×256、256×240、256×256。根据分块结果,按式(12)组装得到新的矩阵和/>维度均为512×512。
步骤6.2:按式(13)定义两个指标矩阵Eindex和E′index,维度均为4×128,分别为
进一步按式(14)组装得到指标矩阵Eindex,维度为8×128,表示为
步骤6.3:根据指标矩阵Eindex,按式(15)对其进行分块,得到到/>以及/>到/>其中/>到/>的维度均为8×4,/>到/>的维度均为8×4。按式(16)对指标矩阵Eindex进行重组,得到新的指标矩阵/>表示为
步骤6.4:将指标矩阵Eindex和拉直,分别记为eindex和/>维度均为1024×1。生成维度为1024×1024的置换矩阵P1,具体方法为:取出eindex和/>的第1个元素1和241,则矩阵P1第241行第1列的元素为1;依次遍历eindex和/>的所有元素,在矩阵P1对应位置填入1,其余元素均为0。
步骤6.5:按式(17)计算C0到Cn的方差ΣD,并按式(18)对矩阵ΣD进行分块,得到矩阵Σ1、Σ2、Σ3、Σ4,维度分别为960×960、960×64、64×960和64×64。
令D=[C0 C1 … Cn-1],计算矩阵M、Q、E、T,按式(19)计算矩阵C0′到Cn′的方差ΣD′。
步骤6.6:按如下步骤生成维度为14400×14400的置换矩阵P2,具体步骤如下:
步骤6.6.1:生成指标矩阵
将拉直,记为/>取出/>的第1个元素0,因为/>则置换矩阵P2的第1行均为0;反之,取出/>的第113个元素1,由于/>则置换矩阵P2第113行第1列元素值为1,其余元素为0。按上述规则遍历/>所有元素,生成置换矩阵P2。按式(21)计算友矩阵Γ的方差矩阵ΣΓ。
步骤7:计算结构模态频率fl、阻尼比ξl和模态振型的方差。
步骤7.1:计算友矩阵Γ的左特征向量ψl,按式(22)计算特征值λl和特征向量相对于友矩阵Γ的灵敏度矩阵/>和/>
步骤7.2:计算再按式(23)计算模态频率fl和阻尼比ξl相对于特征值λl的灵敏度矩阵/>
步骤7.3:通过取矩阵的前8行,得到/>按式(24)计算模态振型/>对于特征向量/>的灵敏度矩阵/>
步骤7.4:按式(25)计算方差矩阵和/>取矩阵/>第一行第一列的元素,即为模态频率fl的方差,取/>第二行第二列的元素,即为阻尼比ξl的方差,/>即为模态振型/>的方差。/>
本发明所公开技术辨识得到的桁架结构模态频率和阻尼比及方差如表1和表2所示,表中“-”表示未辨识出结果,桁架结构模态振型及方差如图4所示。
由表1、表2和图4可以看出,本发明所公开技术可辨识出桁架结构的8阶模态频率、阻尼比及模态振型,并能计算出所有结构参数的方差。计算得到的参数方差与参数相比很小,说明辨识出的结构参数可靠性很高,并且这些结构参数受到测量噪声影响很小。同时,在8阶模态中,与其他阶模态相比,第六阶和第七阶模态频率、阻尼比及模态振型的方差更大,说明第六阶和第七阶模态受到测量噪声的影响更大,辨识结果可靠性较低。在本实施例中,造成该现象的原因是第六和第七阶模态振动幅值很小,属于弱模态。
通过上述分析,可以看出本发明所公开方法可有效辨识出结构参数,并且通过计算得到的方差可对辨识结果的可靠性进行评估。
表1本发明公开技术辨识得到的桁架结构模态频率及方差
表2本发明公开技术辨识得到的桁架结构阻尼比及方差
为突出本发明所公开技术的优势,采用发明公开号CN107561934 A公开的技术辨识出的桁架结构模态频率如表3所示,采用已有的频域分解技术辨识出的桁架结构模态频率如表4所示,表中“-”表示所使用的技术不能辨识出结果。
表3发明公开号CN107561934 A公开的技术辨识得到的桁架结构模态频率
/>
表4已有的频域分解技术辨识得到的桁架结构模态频率
通过对比表1和表3、表4可以看出,本发明所公开方法可解决非白噪声载荷下的结构参数辨识问题,辨识出的结构模态频率精度与真实值更接近,并且可计算结构参数的方差。然而,已有的发明公开号CN 107561934 A公开的技术辨识出的结构参数虽然不包含虚假模态,但辨识出的模态频率精度较低,特别是第八阶模态频率误差为0.26Hz,远大于本发明所公开方法得到的误差0.09Hz,且该已有技术无法计算结构参数的方差;已有的频域分解技术包含虚假模态,应用于结构健康监测、状态评估时会导致误判,此外,该已有技术辨识出的模态频率误差很大,例如第七阶模态频率误差为0.3Hz,而且无法计算结构参数的方差。因此,与已有技术相比,本发明公开的方法具有显著优势。
步骤8:步骤7计算得到的模态频率fl、阻尼比ξl和模态振型的方差可用于评估步骤5中辨识的模态频率fl、阻尼比ξl和模态振型/>的精度,方差越小,说明辨识精度越高,结构振动响应测量噪声越小。进一步,将辨识的结构模态频率fl、阻尼比ξl和模态振型/>及其方差应用于工程结构的健康监测及状态评估。
针对目前工作模态分析方法在实际工程应用时存在的上述问题,本发明公开的基于参数化传递率的结构参数辨识及方差计算方法要解决的问题是,利用两个工况下的相关函数,构建参数化模型辨识结构参数,并计算辨识的结构参数的方差。同时,本发明所公开的方法可处理非白噪声载荷下的结构参数辨识问题,辨识结果中不包含非白噪声载荷引起的虚假模态,并且根据计算得到的方差可对辨识的结构参数精度进行评估,具有重要的工程应用意义。
本发明所公开的方法首先计算结构振动响应的相关函数,随后利用相关函数构建参数化模型,因此区别于发明公开号CN107561934 A中直接利用振动响应构建参数化模型的技术,此外本发明公开的方法可计算结构参数的方差,具有更强的工程实用性。
虽然,上文中已经用一般性说明及具体实施例对本发明作了详尽的描述,但在本发明基础上,可以对之作一些修改或改进,这对本领域技术人员而言是显而易见的。因此,在不偏离本发明精神的基础上所做的这些修改或改进,均属于本发明要求保护的范围。
Claims (10)
1.基于参数化时域传递率的结构参数辨识及方差计算方法,其特征在于:包括以下步骤:
步骤1:采用位移传感器、速度传感器或加速度传感器,测量两个不同工况下工程结构的位移、速度或加速度振动响应信号;
步骤2:确定转换通道,计算结构振动响应信号的相关函数;
步骤3:确定相关函数中的参考信号和非参考信号,在两个工况下构建参数化模型,估计出参数化模型的系数矩阵;
步骤4:计算系数矩阵的方差矩阵;
步骤5:根据步骤3得到的系数矩阵,构造友矩阵,再计算结构模态频率、阻尼比和模态振型;
步骤6:计算友矩阵的方差矩阵;
步骤7:计算结构模态频率、阻尼比和模态振型的方差。
2.根据权利要求1所述的基于参数化时域传递率的结构参数辨识及方差计算方法,其特征在于:步骤1中,采用加速度振动响应信号,两个工况下测量的振动响应信号分别记为x(1)[t]和x(2)[t],其中x(1)[t]和x(2)[t]均是维度为No×1的向量,上标(1)和(2)表示第一个工况和第二个工况,[t]表示第t个时刻,No表示传感器个数并且为大于等于2的整数。
3.根据权利要求1所述的基于参数化时域传递率的结构参数辨识及方差计算方法,其特征在于:步骤2的具体包括:
步骤2.1:根据步骤1得到的振动响应信号x(1)[t],选取x(1)[t]中均方根最大的元素,记为其中z为整数且1≤z≤No,表示振动响应信号向量x(1)[t]的第z个元素,则转换通道为z;
步骤2.2:采用公知的相关函数估计方法,分别计算振动响应x(1)[t]与振动响应x(2)[t]与/>的相关函数,分别记为s(1)[t]和s(2)[t],其中,s(1)[t]和s(2)[t]均为维度是No×1的向量。
4.根据权利要求1所述的基于参数化时域传递率的结构参数辨识及方差计算方法,其特征在于:步骤3的具体包括:
步骤3.1:在No个传感器中挑选出Nr个传感器,其中Nr为整数并且 表示向下取整;将挑选出的Nr个传感器作为参考传感器,对应的相关函数为参考信号,两个工况下的参考信号分别记为/>和/>其余的Nu=No-Nr个传感器为非参考传感器,对应的相关函数为非参考信号,分别记为/>和/>
步骤3.2:根据步骤3.1得到的和/>构建参数化模型如式(1)所示:
式(1)中,上标“(k)”表示第k个工况,和/>表示第k个工况下的系数矩阵,e(k)[t]表示第k个工况下的残差序列,n表示模型阶数,模型阶数n可利用公知的赤池信息准则进行确定;
步骤3.3:式(1)中,所有的系数矩阵和/>由式(2)计算得到
式(2)中,上标“T”和“-1”分别表示矩阵转置和求逆,N为时刻点个数,Θ(k)和Ψ[t]分别为
通过式(2)计算得到Θ(k)后,即可按照式(3)的第一式得到系数矩阵和/> 和的矩阵维度分别为Nu×Nu和Nu×Nr。
5.根据权利要求1所述的基于参数化时域传递率的结构参数辨识及方差计算方法,其特征在于:步骤4的具体包括:
步骤4.1:根据式(2)计算得到的Θ(k),计算两个工况下的残差序列e(k)[t],如式(4)所示
步骤4.2:根据式(4)计算得到的e(k)[t],计算两个工况下残差序列e(k)[t]的方差矩阵,如式(5)所示
步骤4.3:计算两个工况下系数矩阵Θ(k)的方差矩阵,如式(6)所示
6.根据权利要求1所述的基于参数化时域传递率的结构参数辨识及方差计算方法,其特征在于:步骤5的具体包括:
步骤5.1:令为维度Nu×Nu的单位矩阵,即/>其中/>表示维度Nu×Nu的单位矩阵,构造矩阵Ci如式(7)所示
进一步,计算其中上标“+”表示求矩阵伪逆,随后构造友矩阵Γ如式(8)所示
式(8)中,和/>分别表示维度为No×No的单位矩阵和全零矩阵;
步骤5.2:对友矩阵Γ进行特征值分解,记得到的第l个特征值及对应的特征向量分别为λl和则第l阶模态频率fl和阻尼比ξl分别为
式(9)中,Ts为振动响应信号采样的时间间隔,ln(·)和|·|分别表示取自然对数和求绝对值运算;
步骤5.3:根据特征向量计算结构模态振型向量/>首先取特征向量/>的前No个元素,记为/>随后按式(10)计算模态振型向量
式(10)中,为向量/>的第m个元素,Re(·)表示取实部,m的取值是任意的,只需满足即可,优选/>值较大的元素;
由步骤5.2和步骤5.3计算得到的模态频率、阻尼比和模态振型即为最终辨识的结构参数。
7.根据权利要求1所述的基于参数化时域传递率的结构参数辨识及方差计算方法,其特征在于:步骤6的具体包括:
步骤6.1:将步骤4得到的系数矩阵Θ(k)的方差矩阵进行分块,如式(11)所示
式(11)中,的维度为/>的维度为(n+1)NuNr×(n+1)NuNr,/>与/>的维度分别为/>和/>根据分块结果,组装得到新的矩阵/>如式(12)所示
步骤6.2:定义两个指标矩阵Eindex和E′index,如式(13)所示
E″index=(n+1)NuNo+Eindex
进一步组装得到指标矩阵Eindex,如式(14)所示
步骤6.3:根据步骤6.2得到的指标矩阵Eindex,对其进行分块,如式(15)所示
式(15)中,每个矩阵的维度均为2Nu×Nu,每个/>矩阵的维度均为2Nu×Nr,根据分块结果,对指标矩阵Eindex进行重组,得到新的指标矩阵/>如式(16)所示
步骤6.4:将步骤6.2和6.3得到的指标矩阵Eindex和拉直,分别记为eindex和/>生成置换矩阵P1,其维度为2(n+1)NuNo×2(n+1)NuNo;生成置换矩阵P1的具体方法为:取出eindex的第ε个元素,记为eε,取出/>的第ε个元素记为/>矩阵P1第/>行第eε列的元素为1,整数ε从1到2(n+1)NuNo依次递增,矩阵P1其余元素均为0;
步骤6.5:计算C0到Cn的方差ΣD,如式(17)所示
对矩阵ΣD进行分块,如式(18)所示
式(18)中,Σ1矩阵的维度为2nNuNo×2nNuNo,Σ4矩阵的维度为2NuNo×2NuNo;
计算矩阵C′0到C′n的方差ΣD′,如式(19)所示
ΣD′=(MT+Q)Σ4(MT+Q)T+EΣ1ET+(MT+Q)Σ3ET+EΣ2(MT+Q)T (19)
式(19)中,其中D=[C0 C1 … Cn-1];T矩阵满足/>其中vec(·)表示矩阵拉直运算;
步骤6.6:生成维度为(nNo)2×(nNo)2的置换矩阵P2,具体步骤如下:
步骤6.6.1:生成指标矩阵如式(20)所示
式中表示维度为(n-1)No×nNo的零矩阵;将/>拉直,记为/>
取出的第ε个元素,记为/>若/>则置换矩阵P2的第ε行均为0,反之,若/>则置换矩阵P2第ε行第/>列元素值为1;遍历整数ε从1到2(n+1)NuNo依次递增,矩阵P2其余元素均为0;
计算友矩阵Γ的方差矩阵ΣΓ,如式(21)所示
8.根据权利要求1所述的基于参数化时域传递率的结构参数辨识及方差计算方法,其特征在于:步骤7的具体包括:
步骤7.1:计算特征值λl和特征向量相对于友矩阵Γ的灵敏度矩阵/>和/>如式(22)所示
式(22)中,ψl为Γ的第l阶左特征向量,即满足
步骤7.2:计算模态频率fl和阻尼比ξl相对于特征值λl的灵敏度矩阵如式(23)所示
式(23)中,Im(·)表示取虚部;
步骤7.3:计算模态振型对于特征向量/>的灵敏度矩阵,如式(24)所示
式(24)中,通过取式(22)得到的/>矩阵的前No行得到,/>和/>分别表示维度为No′(m-1)和No′(No-m)的零矩阵;
步骤7.4:计算模态频率fl、阻尼比ξl和模态振型的方差,如式(25)所示
式(25)中,取矩阵第一行第一列的元素,即为模态频率fl的方差,取/>第二行第二列的元素,即为阻尼比ξl的方差,/>即为模态振型/>的方差。
9.根据权利要求1所述的基于参数化时域传递率的结构参数辨识及方差计算方法,其特征在于:步骤7计算得到的模态频率fl、阻尼比ξl和模态振型的方差可用于评估步骤5中辨识的模态频率fl、阻尼比ξl和模态振型/>的精度。
10.根据权利要求6所述的基于参数化时域传递率的结构参数辨识及方差计算方法,其特征在于:辨识的结构模态频率fl、阻尼比ξl和模态振型及其方差应用于工程结构的健康监测及状态评估。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202311664658.8A CN117669211A (zh) | 2023-12-06 | 2023-12-06 | 基于参数化时域传递率的结构参数辨识及方差计算方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202311664658.8A CN117669211A (zh) | 2023-12-06 | 2023-12-06 | 基于参数化时域传递率的结构参数辨识及方差计算方法 |
Publications (1)
Publication Number | Publication Date |
---|---|
CN117669211A true CN117669211A (zh) | 2024-03-08 |
Family
ID=90074694
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202311664658.8A Pending CN117669211A (zh) | 2023-12-06 | 2023-12-06 | 基于参数化时域传递率的结构参数辨识及方差计算方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN117669211A (zh) |
Citations (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN107561934A (zh) * | 2017-08-24 | 2018-01-09 | 北京理工大学 | 基于多参考时域传递率的仅输出结构模态参数辨识方法 |
CN112798253A (zh) * | 2021-01-20 | 2021-05-14 | 南京航空航天大学 | 一种考虑非白环境载荷影响的结构模态参数辨识方法 |
US20210223214A1 (en) * | 2019-10-24 | 2021-07-22 | Dalian University Of Technology | Method of complex modal identification for the structure with proportional damping |
CN116911049A (zh) * | 2023-07-28 | 2023-10-20 | 南京航空航天大学 | 单段振动响应数据的结构模态参数不确定性量化方法 |
-
2023
- 2023-12-06 CN CN202311664658.8A patent/CN117669211A/zh active Pending
Patent Citations (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN107561934A (zh) * | 2017-08-24 | 2018-01-09 | 北京理工大学 | 基于多参考时域传递率的仅输出结构模态参数辨识方法 |
US20210223214A1 (en) * | 2019-10-24 | 2021-07-22 | Dalian University Of Technology | Method of complex modal identification for the structure with proportional damping |
CN112798253A (zh) * | 2021-01-20 | 2021-05-14 | 南京航空航天大学 | 一种考虑非白环境载荷影响的结构模态参数辨识方法 |
CN116911049A (zh) * | 2023-07-28 | 2023-10-20 | 南京航空航天大学 | 单段振动响应数据的结构模态参数不确定性量化方法 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
Wang et al. | Event-based state estimation under constrained bit rate: An encoding–decoding approach | |
JP5344823B2 (ja) | センサにより取得されたサンプルのストリームにおける変化を検出する方法 | |
Luong et al. | Observability, redundancy, reliability and integrated design of measurement systems | |
CN111767791A (zh) | 一种基于对抗正则化深度神经网络的到达角估计方法 | |
CN111812980B (zh) | 基于未知输入观测器的离散切换系统的鲁棒故障估计方法 | |
CN113051822A (zh) | 基于图注意力网络与lstm自动编码模型的工业系统异常检测方法 | |
CN108009566B (zh) | 一种时空窗口下的改进型pca损伤检测方法 | |
CN112396088A (zh) | 一种小样本下隐式激励对抗训练的机械故障智能诊断方法 | |
CN111308979B (zh) | 基于多率迟延状态空间模型的u控制系统 | |
CN115186762A (zh) | 一种基于dtw-knn算法的发动机异常检测方法及系统 | |
CN115046766A (zh) | 基于二维灰度图像自适应子空间的小样本轴承故障诊断方法 | |
CN117669211A (zh) | 基于参数化时域传递率的结构参数辨识及方差计算方法 | |
Xiao et al. | Anomalous IoT sensor data detection: An efficient approach enabled by nonlinear frequency-domain graph analysis | |
CN114091600A (zh) | 一种数据驱动的卫星关联故障传播路径辨识方法及系统 | |
CN110662220B (zh) | 基于时空相关性和信息熵的无线传感器网络异常检测方法 | |
CN116362140A (zh) | 基于Transformer和GRU融合的传感器误差补偿方法、系统、装置及介质 | |
CN116992224A (zh) | 一种基于多头注意力机制的时间序列数据重构方法 | |
CN109214268B (zh) | 一种基于集成流形学习的填料塔液泛状态在线监测方法 | |
CN115422603A (zh) | 一种基于模板主成分回归的分组密码硬件安全评估方法 | |
CN115482446A (zh) | 基于深度卷积生成对抗网络的轨迹异常检测方法 | |
Wu et al. | A switching strategy for robust cooperative exploration | |
Kawada et al. | An application of a generalized inverse regression estimator to Taguchi's T-Method | |
CN109274107B (zh) | 一种计及奇异值的低频振荡信号参数辨识方法 | |
CN113792776A (zh) | 面向网络安全异常检测中深度学习模型的解释方法 | |
CN112883641B (zh) | 一种基于优化elm算法的高大建筑倾斜预警方法 |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination |