Python开源模块：CESSIPy - 用于土木工程随机系统识别的Python工具

≥软件X 18（2022）101091原始软件出版物CESSIPy：一个用于土木工程随机系统识别的Python开源模块Matheus Roman Carini，Marcelo Maia Rocha巴西阿雷格里港南里奥格兰德联邦大学ar t i cl e i nf o文章历史记录：收到2021年收到修订版，2022年3月16日接受，2022年保留字：系统识别工作模态分析结构动力学增强频域分解挑峰随机子空间识别a b st ra ct本文介绍了CESSIPy，这是一个Python开源模块，用于从仅输出的测量中识别振动结构的模态特性。识别的属性是自然频率，阻尼比和模态振型。假设系统输入为宽带平稳高斯过程。现有五种系统辨识方法：（1）基本频域（BFD），（2）增强频域分解（EFDD），（3）协方差驱动随机子空间辨识（SSI COV），（4）数据驱动随机子空间辨识（SSI DATA），（5）辅助变量（IV）。所有功能都通过数值模拟进行了验证。该软件有完整的文档，包括用户指南、代码内描述和脚本示例，以演示代码的使用。©2022作者（S）。由爱思唯尔公司出版这是CC BY许可下的开放获取文章（http://creativecommons.org/licenses/by/4.0/）中找到。代码元数据当前代码版本v1.1此代码版本使用的代码/存储库的永久链接https://github.com/ElsevierSoftwareX/SOFTX-D-21-00072Code Ocean compute capsule法律代码许可证MIT许可证使用的代码版本控制系统使用Python 3的软件代码语言、工具和服务（模块Matplotlib、MRPy、NumPy和SciPy）编译要求、操作环境依赖性Python 3.6如果可用，链接到开发人员文档/手册http://github.com/MatheusCarini/CESSIPy问题支持电子邮件matheuscarini@gmail.com1. 动机和意义系统辨识的目的是从观测数据中拟合系统动态行为的数学模型中所需的参数[1]。在土木工程中，桥梁和塔架等结构是建模的结构系统，相应的模型参数是特征频率、阻尼比和模态振型[2]。对于大型结构，即使不是不可能，也很难进行强迫振动试验。然后，可以在正常操作条件下应用的识别技术环境测试不需要关闭流量或中断正常运行，但分析收集的数据可能具有挑战性[3]。*通讯作者。电子邮件地址：matheuscarini@gmail.com（Matheus Roman Carini），mmrocha@ufrgs.br（Marcelo Maia Rocha）。https://doi.org/10.1016/j.softx.2022.101091在这些技术中，假设结构由不可测量的输入力激励，并且仅输出测量值（例如，加速度）。因此，随机识别问题包括仅从给定的输出数据计算结构系统矩阵[4]。载荷通常假定为宽带平稳分布和高斯分布。一个重要的应用是数值模型的验证，这通常是通过有限元方法创建的一旦结构建成，就记录环境振动测量值。接下来，这些测试结果用于验证数值模型或校准其参数[5]。系统辨识技术也应用于结构健康监测。土木工程结构中损伤的累积会引起结构动力特性的变化。通过比较结构在整个寿命周期内的模态参数，可以实现结构损伤的检测。自然频率的变化是最2352-7110/©2022作者。 由Elsevier B.V.出版。这是一篇开放获取的文章，使用CC BY许可证（http://creativecommons.org/licenses/by/4.0/）。可在ScienceDirect上获得目录列表SoftwareX期刊主页：www.elsevier.com/locate/softxMatheus Roman Carini和Marcelo Maia Rocha软件X 18（2022）1010912R⃗R[客户端][客户端]⏐2图1.一、CE S S I P y 模 块 的 示意图。使用损坏指示器[6]，尽管根据作者有几个商业软件包用于随机系统识别，如ARTeMIS [7]。然而，这些软件在学术环境中有一些缺点。例如，不可能访问算法的实现。因此，大多数研究人员开发自己的代码。2.1.1. BFD方法经典的BFD方法-也称为峰值拾取（PP）技术-在许多书籍中描述， Felber [10]开发了一种实用而系统的实施方法。本征频率通过选择输出谱密度的峰值来确定在谐振频率fr处，谱密度矩阵可以近似为：尽管Python应用程序稳步增长，但作者并不知道存在开源代码，Gy（fr）αrr H（二）本文提出的识别方法。Python是一种免费、开源、面向对象的编程语言。它可以在大多数操作系统上运行，并且在科学计算中变得非常流行。CESSIPy （ Civil Engineer Stochastic System Identification forPython）是一个Python模块，用于从加速度测量中估计本征频率、模态形状和阻尼比其中，αr是由第r个自然频率确定的因子，阻尼比和模态振型，而εr是第r阶模态振型估计。本文提出了两种阻尼比的估计方法半功率带宽方法根据峰值频率fr和两个半功率频率fa和fb估计阻尼：的部分。 图 1给出了CESSIPy的示意图。的=f2−f2（三）模块已通过数值模拟验证，实验数据算例及模态分析阿吉布一4f2Z24桥的详细信息见第4节。CESSIPy GitHub存储库上提供了用户指南和脚本示例。所有函数都有代码内注释，描述函数本身以及输入和输出参数。另一种选择是使用曲线拟合程序将峰值附近的响应自动谱调整为由具有恒定谱值的输入激励的单自由度（SDOF）系统的理论fRR2. 系统辨识方法Gy（f）=C1<$1−（f/f）2+j2<$（f/f）<$+C2（4）CESSIPy上提供了五种方法，本节简要介绍了这些方法。在软件用户指南和相关文献中可找到更深入的文本（参见[2，3，8]）。2.1. 频域基本频域（BFD）和增强频域分解（EFDD）方法假设输入是白噪声随机过程。两种算法都需要响应谱密度矩阵Gy（f），其估计为：Gy（f）=EY（f）YH（f）（1）其中Y（f）是频域结构响应。 superindex（.）H表示复共轭转置，E. 表示期望值。其中C1是由系统刚度确定的系数，输入谱，fr为无阻尼固有频率，ω为阻尼比，j为虚部，C2为考虑其它模态影响的因子。2.1.2. EFDD方法EFDD方法由Brincker等人开发[11，12]。 第一步是计算谱矩阵的奇异值分解（SVD）：Gy（f）=U（f）<$（f）UH（f）（5）其中U（f）是包含降序奇异值的对角矩阵，U（f）是包含奇异向量的复酉矩阵。 该矩阵被分解成一组自谱密度函数，每个自谱密度函数对应于一个单自由度系统。在谱峰附近，奇异向量是模态形状的估计。如果只有一种模式是相关的2Matheus Roman Carini和Marcelo Maia Rocha软件X 18（2022）10109131=compM微米p1微米⃗我--⃗⃗=I微米DMDMMM0=考虑到谐振频率fr，频谱矩阵可以近似为考虑最大奇异值的秩一矩阵：Gy（f）<$=σ1（fr）u<$1（fr）u<$H（fr）（6）第一奇异向量u1（fr）是模态振型的估计。峰值附近的功率谱密度函数2.2.2. IV法IV方法评估自回归移动平均（ARMA）模型。然而，仅以线性方式识别AR参数，避免了实现了由Peeters [2]提出的参考传感器的算法导出的公式从p阶的逆时ARMA模型开始：通过逆快速转换器转换回时域。y参考+αby参考+···+αby参考=eb+γbeb+···+γbeb（十六）傅立叶变换（FFT）。 自由衰减时域函数k1k+1pk+pk1k+1p k+p因此获得，这是自相关函数的其中，y=ref是输出参考传感器矢量，而y=b是a单自由度系统其固有频率和阻尼比分别为kk通过曲线拟合自由衰减函数估计：白噪声矢量左侧称为自回归（AR）部分，右侧称为移动平均（MA）部分。矩阵αb是AR矩阵参数，y（t）=ypexp（−ηt） cos（ωdt）（7）γbl其中yp为振幅，η ωnn为振幅，ωn和ωd分别为无阻尼和有阻尼的固有圆频率。2.2. 时域l是MA矩阵参数。后向AR系数由输出协方差通过最小二乘计算。这些系数聚集在一个矩阵中，称为状态矩阵，其形式为Ab，p。矩阵包含 伴随矩阵的特征值和特征向量是分别计算为Λd和Λcomp固有频率和阻尼比由µi特征值计算得出在第二步中确定模态形状。 减少随机模态参与矩阵Gref由协方差驱动（SSI COV）和数据驱动（SSI DATA）mAb，p在Python中，注释-在CESSIPy上实现了随机子空间辨识，的特征向量矩阵的前r行：comp以及工具变量（IV）方法。每种方法的主要步骤如下。Gref=（λcomp[r，：]）T（17）其中r是参考传感器的数量。反前-2.2.1. SSI COV和SSI DATA方法趋向模态随机能控性矩阵评价：SSI技术识别离散时间随机状态空间从仅输出的数据建模该模型由rref=[Λp−1GrefΛp−2Gref. . .ΛdGrefG参考]（十八）以下等式：最后，估计模态振型：xk+1=AV=[Rrefrefp−1. . .R参考]（Γref）−1（19）yk=Cxk+zk（9）其中，xk是在时间k（k）处的维度n的状态向量、1、. . . ，N1，其中N是时间样本的数量），wk是由于干扰和建模不准确性引起的过程噪声，y k是输出向量，z k是由于传感器不准确性引起的测量噪声，A是状态矩阵，C是观测矩阵。通过输出协方差矩阵（SSI COV）或将未来输出的行空间投影到过去输出的行空间（SSI DATA）来估计模型矩阵，并计算特征频率fi、阻尼比i和振型i发件人：A=Λd−1（10）Λd= diag（μi），i = 1，. . . ，n（11）lnµi其中，Rref是在所有输出和参考传感器之间计算的简化输出协方差矩阵。3. 软件描述CESSIPy 在 Python 3 中 实 现 。 它 需 要 MRPy [14] ， 可 以 从GitHub 免 费 获 得 ， 以 及 流 行 的 库 ： NumPy ， Matplotlib 和SciPy。该模块可以从GitHub存储库免费下载它被设计为在脚本模式下或在一个笔记本电脑中运行。3.1. 软件构架CESSIPy模块包含22个功能，如图所示。二、它们被分为三组：频域相关函数、时域相关函数和其他函数。主要功能如下。λi=（十二）Δt3.1.1. 时域相关函数F|λi|2π实数（λi）（十三）标识开始调用rearrange_data函数，该函数排列定位引用的时间数据在输出阵列的第一行中的传感器的迭代器RMatheus Roman Carini和Marcelo Maia Rocha软件X 18（2022）1010914||i=（十四）|λi|SSI_COV_iterator、SSI_DATA_iterator、IV_iterator和V=C=[1. . .你好. .（15）Fast_SSI是返回不同阶模型的本征频率、阻尼比和模态形状的包装函数。其中，λ i和Λd是A的特征向量和特征值矩阵，μi和λi分别是离散和连续特征值，λt是采样时间，λ i是连续特征值。 表示复数模量。Peeters [2] 提出 的基于 参考的 算法在 CESSIPy上实 现。 由Döhler和Mevel [13]导出的快速SSI算法2也可用。它们调用子程序来执行一些后台计算。例如，SSI_COV、SSI_DATA和IV子例程根据给定的模型阶次估计模态模型。迭代器的结果作为stabilization_diagram函数的参数。 顾名思 义 ， 它 绘 制 稳 定 图 并 返 回 稳 定 极 点 。 最后 ， 函 数stable_modes将接近的极点聚集到单个模式中。该过程是自动化的，在第4节中进行了说明。Matheus Roman Carini和Marcelo Maia Rocha软件X 18（2022）10109153.1.2. 频域相关函数频域识别程序通过调用SDM函数开始估计频谱密度矩阵。如果用户需要，它会生成频谱矩阵图（图1显示了一个示例）。谱密度函数是BFD和EFDD算法的参数。BFD和EFDD功能在默认交互模式下运行时都例如，通过在ANPSD_from_SDM函数生成的平均归一化功率谱密度图中单击来选择固有频率。还有一个批处理模式，允许自动识别。在这种模式下，用户应该在调用方法函数之前分配一些属性（在用户指南中有描述）3.1.3. 其他功能一些附加功能不适合在前两个章节。MAC功能计算模态保证标准[15]。相干函数有一个不言自明的名字，它计算并绘制信号之间的相干函数。plot_1dshapes和plot_3das1d函数分别绘制梁状结构和建筑物状结构的模态模型。3.2. 软件功能时域相关功能向用户提供以下计算输出协方差、Toeplitz矩阵和投影矩阵。通过SSI COV、SSI DATA和IV算法估计本征频率、阻尼比和模态振型。• 绘制稳定图。频域相关功能提供以下功能：计算并绘制功率谱密度矩阵、平均归一化功率谱密度和相干函数。通过BFD和EFDD算法估计特征频率、阻尼比和模态振型。4. 应用4.1. 说明性示例源代码1说明了CESSIPy这里只演示时域方法，因为它们是自动化的，因此不需要用户交互。更多示例，包括频域算法，可以在GitHub存储库中找到。模拟结构为图1所示的三层平面框架。1.一、激振力是随机的（ 由numpy.random.randn（） 方 法 产 生 的 零 均 值 高 斯 白 噪声）。只考虑了三个自由度，分别对应于每层楼板的水平位移。首先，对振动实验进行了仿真。特征频率和模态振型由质量和刚度矩阵计算。所有的模态阻尼比，α，设置为等于1.5%。施加随机力，并使用MRPy类计算结构响应。然后，从加速度响应向量yk估计模态模型。请注意，每种方法只需要四个脚本行来估计和绘制模态参数。迭代器返回的三个数组是已识别的自然频率、阻尼比和模态形状（例如FNCR、ZTCR、用于SSI COV方法的VVCR阵列）。该脚本生成四个独立的模态形状图，它们被组装在一起，如图3所示。实际和识别的模态振型几乎相同。源代码1：用于模拟振动实验和演示时域系统识别算法的代码。值得注意的是，由于力的随机性，每次模拟都会导致不同的结果。对100个不同的负载模拟重复源代码1。结果综合在表1和表2中，表1和表2分别给出了100次仿真的平均本征频率和阻尼比估计。所有的方法都能产生很好的本征频率····Matheus Roman Carini和Marcelo Maia Rocha软件X 18（2022）1010916表1图三. 实模态和识别模态形状.图二. CESSIPy模块的功能。两个中央桥墩被夹在主梁上，而两个桥台通过混凝土铰链连接到主梁上[2，16]。测试数据可供科学界使用，并在许多研究中用作基准。在CESSIPy上实现的时域技术应用于Z24桥的振动数据进行验证。本文所用的环境试验数据来源于渐进损伤试验8。在这种情况下，系统响应是在33个加速度计的9个单独记录（设置）中测量的以100Hz的速率对数据进行采样，设置的测量时间为10 min 55.36 s，对应于65，536个样本。分别对每种安装进行模态识别。表3和表4分别列出了平均本征频率和阻尼比估计值，以及Brincker等人[17]、Luscher等人[18]和Peeters等人[19]获得的结果。在0-20 Hz的频率范围内识别出七种模式。对于大多数模态，频率和阻尼估计都与其他作者的结果一致。5. 影响来自100次模拟的实际和平均本征频率估计（Hz）。模式房SSI COVSSI数据IV月15,465,465,465,46第215,2915,2915,2915,29第322,1022,1022,1022,10表2实际和平均阻尼比估计从100模拟。模式房SSI COVSSI数据IV月11.50%1.51%1.49%1.50%第21.50%1.47%1.43%1.48%第31.50%1.47%1.40%1.47%估计。SSI的阻尼比估计值COV和IV非常接近真实值。SSI DATA方法提供的阻尼估计值略低。4.2. Z24桥在欧洲项目SIMCES的框架内，对预应力混凝土桥梁进行了广泛的仪器和测试，目的是检测、定位和量化人为施加的损伤。Z24是一座直桥，略微倾斜。它位于瑞士伯尔尼和苏黎世之间。这是一座后张混凝土箱梁桥，有两条车道，三跨：两个14米的边跨和30米的主跨土木工程结构动力测试的应用是十分广泛的。系统识别提供了可靠的数据，以实验验证结构的数值模型，已在设计阶段。此外，阻尼系统的最终设计需要准确识别结构的动态特性。除此之外，模态特性的识别在发生修改、修复或强化干预时起着重要作用[20]。然而，缺乏用于计算实现、改编和使用的免费和开源代码。CESSIPy包含五个著名的随机系统识别方法。该模块可以显着帮助学生，研究人员和专业人士以简单可靠的方式从仅输出数据中识别模态模型。例如，研究人员可以在自己的软件中改编和实现代码。工程师可以很容易地从加速度数据估计模态模型，而无需深入了解方法最后，CESSIPy作者将其免费提供下载，以增加模块的可见性和使用率。考虑到Python令人难以置信的增长和对系统识别越来越感兴趣，作者预计CESSIPy将对结构工程师非常有用。Matheus Roman Carini和Marcelo Maia Rocha软件X 18（2022）1010917表3本征频率估计值的比较（Hz）。模式SSI COVSSI数据IVBrincker等人[17]Luscher等人[18个国家]Peeters等人[19个]13,863,863,863,863,853,85924,904,914,904,894,894,9039,759,759,769,789,829,77410,2910,2910,2910,3110,3010,31512,4012,4212,4012,5012,7012,5613,2913,3313,2613,10–13,19719,1619,0619,04––19,0表4阻尼比估计值的比较（%）。模式SSI COVSSI数据IVBrincker等人[17]Luscher等人[18个国家]Peeters等人[19个]10,840,950,880,790,800,921,261,241,421,611,521,431,361,341,391,261,831,341,461,581,501,361,741,452,902,862,851,366,052,564,143,673,401,65–3,072,192,282,14––2,06. 结论本文介绍了CESSIPy，这是一个开源的Python模块，可以从仅输出的测量值中估计模态模型。三种时域方法可用：SSI COV，SSI DATA和IV;两种频域方法：BFD和EFDD。此外，还创建了用户指南、代码内描述和脚本示例，以演示应如何使用和实现代码。CESSIPy是土木工程结构随机识别的有力工具该模块可用于教育，研究和专业目的。竞合利益作者声明，他们没有已知的竞争性财务利益或个人关系，可能会影响本文报告的工作致谢作者要感谢KU Leuven Struc- tural Mechanics Section分享Z24桥梁振动数据，这些数据记录在欧洲项目SIMCES的框架这项工作得到了巴西国家科学技术发展委员会（CNPq）和巴西高等教育人员改进协调会（CAPES）的支持。引用[1] 永湖系统识别：用户的理论。第二版 UpperSaddleRiver，USA：PrenticeHall，1999.[2] 彼得斯湾土木工程中的系统识别与损伤检测。（博士）Leuven：KatholiekeUniversiteit Leuven; 2000.[3] Brincker R，Ventura CE. 操作模态分析导论。奇切斯特：约翰威利父子;2015年，http://dx.doi.org/10.1002/9781118535141。[4] 范奥弗施，德穆尔.线性系统的子空间辨识：理论、实作与应用。Boston：Kluwer; 1996，http://dx.doi. org/10.1007/978-1-4613-0465-4。[5] Juang J-N，Pappa RS.模态参数识别和模型降阶的特征系统实现算法。JGuid Control Dyn 1985;8（5）：620-7. http://dx.doi.org/10.2514/3.20031网站。[6]Andersen P.使用向量ARMA模型识别土木工程结构。（博士）论文），奥尔堡：奥尔堡大学; 1997。[7] 结构振动解决方案。ARTEMIS。2020年，https://svibs.com/。[2020年10月17日[8] Rodrigues J. Identificação modal estocástica：metodos de análise e apli-caçées em estruturas de engenharia civil.（博士）论文），波尔图，葡萄牙：波尔图大学; 2004年，https://hdl.handle.net/10216/11964。[9]Bendat JS，Piersol AG.相关与谱分析的工程应用。第三版New York：JohnWiley and Sons.[10]Felber AJ.混合式桥梁评估系统的开发。 （博士） 论文），温哥华：不列颠哥伦比亚大学; 1993年，http://dx.doi.org/10.14288/1.0050403。[11]Brincker R，Zhang L，Andersen P.基于频域分解的环境响应模态识别。第18届国际模态分析会议论文集：2000年2月7-10日;美国圣安东尼奥。1999年，p.625-30[12]Brincker R，Ventura CE，Andersen P.通过频域分解进行阻尼估计。第19届国际模态分析会议论文集：2001年2月5-8日;美国基西米。1999，p. 698-703[13] Döhler M，Mevel L. 基于子空间的系统辨识中系统矩阵的快速多阶计算。控制工程实践2012;20（9）：882-94。http://dx.doi.org/10.1016/j.conengprac.2012.05.005网站。[14]RochaMM.Multivariaterandomprocesseswithpython.2020 年 ，github.com/mmaiarocha/MRPy。[15]Allemang RJ，Brown DL.模态向量分析的相关系数。第一届国际模态分析会议论文集：1982年 11 月 8日至10日;美 国 奥兰多。1982年，p.110-6[16] Maeck J，De Roeck G. Z24基准测试的描述。机械系统信号学报2003;17（1）：127-31. http://dx.doi.org/10.1006/mssp.2002.1548网站。[17]Brincker R，Andersen P，Cantieni R. Z24公路大桥的识别与一级损伤检测Exp Tech 2001;25（6）：51网址：//dx.doi.org/10.1111/j.1747-1567.2001.tb00047.x网站。[18]Luscher DJ，Brownjohn JMW，Sohn H，Farrar CR. Z24桥实测数据的模态参数提取。第19届国际模态分析会议论文集：2001年2月5-8日;美国基西米2001年，https：//digital.library.unt.edu/ark:/67531/metadc715804/网站。[19]放大图片作者：Peter B，Maeck J，De Roeck G.土木工程中的激励源与动态系统识别。在：欧洲COSTF3系统识别和结构健康监测会议的会议记录：2000年6月;西 班 牙 马德里。1999年，p.341比50[20]Cunha Jiang，Caetano E，Magalheses F，Moutinho C.桥梁的动态识别和连续动态监测：桥梁生命周期中的不同应用。结构工程2018;14（4）：445-67。http://dx.doi的网站。org/10.1080/15732479.2017.1406959.