基于广义子空间溯踪的时变结构模态参数辨识方法与流程
本发明具体涉及一种基于广义子空间溯踪的时变结构模态参数辨识方法。
背景技术：
水利工程，如大坝的强震观测是利用传感器测量结构在地震激励下实际的震动反应。基于强震观测进行混凝土拱坝模态参数识别，可为结构的抗震分析、健康诊断和震后结构损伤评估等提供基础。基于强震观测来进行混凝土大坝的模态参数识别属于结构的运行模态分析(oma)问题研究的范畴。模态参数，如固有频率、阻尼比、模态质量、模态刚度和模态振型都有明显的物理意义，最能反映大坝结构真实的动力特性。
获取模态参数的途径有两种“理论计算和试验模态分析。理论计算一般是先将实际结构离散成有限自由度系统，而后利用力学和数学方法建立系统振动模型，再运用线性代数、矩阵分析和状态空间等数学方法求解系统状态参数。通过这种途径获取参数方便易行。但由于建模过程必不可少的假设和简化，边界条件及内部结构联结部位又难于实际情况很好吻合，因此可能存在较大误差。试验模态分析是从结构的某些测点上采集的动态数据出发。利用一些信号分析手段和模态参数识别理论来识别模态参数。由于建立在实测的基础之上，显然通过这种途径获得的模态参数能更好的反映结构的实际动力特性。模态参数识别的主要任务就是从测试所得的数据中确定振动系统的模态参数。
基于强震观测的混凝土坝模态识别的研究在国内外已取得相当多的成果。如，h.l.kou,f.jin,j.yng,etal.,modalparameteridentificationofertanarchdamfromstrongearthquakerecords,journalofhydroelectricengineering28(5)(2009)51-56，采用带输入的自回归arx(auto-regressiveexogenous，arx)模型，分别根据二滩拱坝和水口重力坝的强震观测识别了结构的模态参数；l.f.zhang,l.g.xing,strongearthquakeanalysisonlongyangxiagravityarchdam,waterpower12(1998)14-17，根据龙羊峡的强震观测数据识别了结构的模态参数，并和原型激振试验的结果进行了对比；loh等采用随机子空间识别(stochasticsubspaceidentification，ssi)方法，基于强震观测资料对翡翠拱坝的模态参数进行了识别，并研究了库水位对模态参数的影响和非均匀输入对拱坝动力响应的影响等问题。n.okuma,y.etou,k.kanazawa,etal.,dynamicpropertyofalargearchdamafterforty-fouryearsofcompletion,in:proceedingsofthefourteenthworldconferenceonearthquakeengineering,beijing,china,october12-17,2008，根据大坝的振动测试数据，分别采用频域法研究了mauvoisin拱坝和hitotsuse拱坝的系统识别问题。
在本课题组之前的研究中，混凝土重力坝的模态参数问题被看作一个盲源分离(bbs)问题来求解，l.cheng,d.j.zheng,theidentificationofadam'smodalparametersunderrandomsupportexcitationbasedonthehankelmatrixjointapproximatediagonalizationtechnique,mechanicalsystemsandsignalprocessing42(1-2)(2014)42-57。上述研究中采用的模态识别的方法都是基于结构线性振动理论的，即隐含地假定“大坝-地基-水库”系统在一次振动测量过程中保持不变。在强震激励下，由于材料、边界条件、阻尼影响和“大坝-地基-水库”系统之间的相互作用以及结构损伤等因素的影响，混凝土坝可能会出现时变特性。在地震期间坝体混凝土和基岩材料的动力参数可能会发生变化，坝体的接缝可能会随着地震激励振幅的变化而开合。“大坝—水库”和“大坝—地基”之间的相互作用也可能将一些时变因素引入到混凝土坝的动力特性中。除此之外，大坝结构的损伤通常会引起结构刚度、阻尼随时间变化的非线性力学特征，导致结构模态参数出现时变特征。因此，使用定常的结构模态识别方法得到的恒定不变的模态参数，难以实时反映出结构在强震激励下时变性质。
因此，混凝土坝的时变模态识别方法是一个具有重要意义的研究方向。
对于时变结构来说，经典的定常系统模态识别方法的一些假设已不再适用。为此，liu等人对于时变系统基于时间“冻结”的思想，提出了“伪模态”的概念。目前，针对结构时变的模态参数识别问题已发展出了三类方法：时间序列模型法、递归随机子空间法(recursivestochasticsubspaceidentification，rssi)和时频分解法三种方法。公茂盛等采用自适应递归算法(adaptiveforgettingthroughmultiplemodels，afmm)追踪某高层楼房在强震激励下模态参数的变化，并根据模态识别结果进行了结构健康评估。afmm算法可以精确的判断出系统发生时变和跳出时变的时刻，根据这些时刻，整个振动响应记录可以划分为不同的阶段。对于混凝土坝，r.tarinejad,modalidentificationofanarchdamduringvariousearthquakes,internationalconferenceondamsandhydropower,2012通过把pacoima拱坝整个强震记录划分成若干个片段，然后在每一段进行线性定常模态识别，研究了该坝的时变特性。该方法虽然引入时变模态参数识别，但是因识别方法的限制，其对频率、振型的识别精度不足，且方法的鲁棒性有待提高。
“大坝-地基-水库”系统一般具有结构体积较大、运行环境复杂、结构系统自由度高和模态密集等特点。此外，由于受外部干扰较大，混凝土坝的振动测量数据通常具有较低的信噪比。为了追踪强震激励下混凝土坝的时变特性，实现在线的结构健康诊断和损伤评估，研究基于强震观测的混凝土坝时变模态识别方法，具有重要意义。
结构的时变模态识别问题可表述为一个子空间溯踪问题。鉴于此，本发明基于强震观测和rssi方法来研究混凝土坝等土木、水利工程的时变模态参数问题，将子空间溯踪的优势方法—广义子空间溯踪(gyast)引入到rssi算法中，提出了一种新的时变模态参数识别方法。
技术实现要素：
为了解决现有技术中对于结构体积较大、运行环境复杂、模态密集，受外部干扰较大的时变结构，使用定常的结构模态识别方法得到的恒定不变的模态参数，难以实时反映出结构在强震激励下时变性质，其振动测量数据通常具有较低的信噪比，且对频率、振型的识别精度不足，方法的鲁棒性不高等缺陷，本发明提出了一种新的时变模态参数辨识方法。本发明的时变结构模态参数辨识方法可以减小观测噪音和随机激励带来的干扰，鲁棒性好，并且对频率、振型的识别精度较好，此外本方法辨识模态参数的计算效率高，可以更高效的为结构的抗震分析、健康诊断和震后结构损伤评估提供基础。
本发明要解决的技术问题通过以下技术方案实现：
基于广义子空间溯踪的时变结构模态参数辨识方法，其包括以下步骤：
(一)对具有时变特性的结构在正常运行情况下测量得到的结构的振动响应，根据初始一定长度的振动响应监测数据构造多维的时延向量，计算该时延向量的相关函数矩阵，并采用特征值分解来对广义子空间溯踪算法进行初始化；
(二)引入遗忘因子，根据新的振动响应监测数据来对多维时延向量的相关函数矩阵进行在线更新，并采用广义子空间溯踪算法来对正交矩阵π(k)进行在线更新；
(三)最后根据子空间溯踪算法获取的更新的正交矩阵π(k)，计算时变结构任意时刻的离散状态空间矩阵和观测矩阵，并对离散状态空间矩阵进行特征值分析，提取时变结构每一观测时刻的模态参数：模态频率、阻尼比和振型。
进一步地，本发明所述的基于广义子空间溯踪的时变结构模态参数辨识方法中，所述步骤(一)具体为：
步骤(1)，根据被辨识时变结构的工作状态和时变结构的主要频率范围，设定辨识所需的采样频率，并对辨识结构的结构动力学振动响应信号进行采集，得到l个通道的离散的振动观测数据组成的l维向量y0(k)，(k＝1,2,….,n)；k是采样点的序号，n是观测样本总数；
步骤(2)，对时变模态参数识别程序中的参数，系统阶次n，时延p、遗忘因子β和初始数据长度n1进行初始化；
步骤(3)，将步骤(1)采集的振动响应信号，用希尔伯特变换获得相位偏移90度的虚拟响应y90(k)，构造以下的数据向量
并建立由该时延向量组成的矩阵
式(2)中，p为时延，满足p×2l≥n，l表示时延向量的维数，上标t表示转置，n表示系统阶次，
令k＝n1，根据式(2)得到时延向量组成的矩阵y(k)，然后，计算其相关函数矩阵h0(k)，得到的hankel矩阵h0(k)是一个方阵：
对相关函数矩阵h0(k)进行以下的特征值分解；
式(4)中，矩阵q1(k)包含了h0(k)的前n个特征向量，矩阵q1(k)是振动信号的主子空间，代表不同阶模态的贡献，q2(k)是剩余特征值对应的特征向量，δ1(k)和δ2(k)分别是由前n个主特征值和剩余特征值组成的对角矩阵，n表示系统阶次，上标t表示转置；
对正交矩阵π(k)、加速模态度响应的相关函数矩阵
π(k)＝q1(k)；(5)
σn(k)2＝h0(k)的最小特征值(7)。
进一步地，本发明所述的基于广义子空间溯踪的时变结构模态参数辨识方法中，所述步骤(二)具体为：
步骤(4)，令k＝k+1，在tk时刻，获得新的观测数据y0(k)后，通过引入遗忘因子β对hankel矩阵h0(k)进行在线更新，hankel矩阵h0(k)在t＝tk时刻可以更新为
h0(k)＝βh0(k-1)+y(k)y(k)t(8)；
步骤(5)，为了计算加速模态度响应
对于广义子空间溯踪算法，首先引入以下式(9)的近似：
式(9)中，上标h表示复共轭转置，
这时，增广的(n+1)维列空间span(v(k))的标准正交基π(k)可表示为：
π(k)＝[π(k-1),u(k)](10)
式(10)中，u(k)是y⊥(k)的单位范数矢量，
定义矩阵
式(11)中，
式(12)中，q是向量qn(k)的前n个元素组成的向量，rn(k)是向量qn(k)的最后一个元素；
这时可以按照式(13)～(22)计算以下的参数b(k)、s(k)、c(k)、λm(k)和σn(k)2，向量e′(k)和q′(k)，以及矩阵π(k)、c2(k)，以实现对正交矩阵π(k)的在线更新：
b(k)＝1/(1+|rn(k)|)(13)
s(k)＝rn(k)/(|rn(k)|σ(k))(14)
π(k)＝π(k-1)-e′(k)qh(16)
σn(k+1)2＝min{λm(k),σn(k)2}(22)，
参数b(k)、s(k)、c(k)、λm(k)和σn(k)2，向量e′(k)和q′(k)，以及矩阵π(k)、c2(k)均属于计算中的中间计算结果，没有明确的物理含义。
进一步地，本发明所述的基于广义子空间溯踪的时变结构模态参数辨识方法中，所述步骤(三)具体为：
步骤(6)，按照上述步骤(5)的过程对正交矩阵π(k)进行更新后，令矩阵ο＝π(k)，在第k个时间间隔中，离散状态空间矩阵ak和观测矩阵gk使用下面的公式(23)和(24)来进行估计：
gk＝o(1:2l,:)(23)
式(23)和(24)中ο由矩阵o除掉最后的2l行元素后剩余的元素组成，矩阵
得到离散状态空间矩阵ak之后，固有频率和阻尼比使用其特征值进行提取，时变结构的固有模态频率fi从下式获得：
式(25)中，n0是模态阶数，对复模态系统，系统阶数n＝2n0，λi(i＝1,2,…,n0)是ak的特征值，re{·}和im{·}分别表示λi的实部和虚部；
时变结构的阻尼比ξi(i＝1,2,…,n0)使用下式(26)获得：
这时l个测量自由度对应的振型的分量φl，可以表达为
φl＝gkφ(27)
式(27)中，矩阵φ＝[φ1,φ2,...,φn]由ak的n个特征向量组成，φn是ak的第n个特征向量；
步骤(7)，继续运行上述计算步骤(4)～步骤(6)，直至所有的振动观测数据都已分析完成。
进一步地，本发明所述的基于广义子空间溯踪的时变结构模态参数辨识方法中，步骤(1)中，设定辨识所需的采样频率时，采样频率设定满足香农采样定理；
步骤(2)中，对时变模态参数识别程序中的参数，系统阶次n，遗忘因子β和初始数据长度n1进行初始化时，系统阶次n通过绘制奇异值谱图来确定，选定奇异值谱图中第一个转折处对应的奇异值数为系统阶数n；
遗忘因子β的取值根据时变结构的具体工况确定，β取值介于0.95到1.0之间，初始数据长度n1取值为大于或等于400，n1越大模态参数识别越准确，同时计算时间会延长；
时延p满足p×2l≥n，其中，n为系统阶次，l为时延向量的维数。
进一步地，本发明所述的基于广义子空间溯踪的时变结构模态参数辨识方法中，步骤(2)中，对时变模态参数识别程序中的遗忘因子β进行初始化时，遗忘因子β的取值为0.995，初始数据长度n1取值为400。
本发明的有益效果：
本发明的时变结构模态参数辨识方法中，在时变结构模态参数辨识过程中引入“时间冻结”的假设。“时间冻结”假设是将一个线性时变结构用一系列分段时不变的线性定常系统来近似模拟。也就是说，结构的系统参数会在一些离散的时刻突然改变，但在这些时刻之间的小的区间中，保持不变。如果时间间隔划分的足够小，这种分段近似的精度是能够满足工程应用要求的。本发明基于这种思想，采用广义子空间溯踪和递归随机子空间法(gyast-rssi)计算得到离散状态空间矩阵，进而得到每个瞬时的模态参数：固有频率、阻尼比和模态振型。本发明从参数化的时频域的角度出发，给出一种基于时频分析的参数化时域的时变结构模态参数辨识方法，物理意义清晰；适用于土木、水利工程应用领域的时变结构的模态参数辨识，并且具有使用简单方便、计算精度高、鲁棒性强和计算效率高等特点。
本发明的时变结构模态参数辨识方法用于水利工程时变结构，如混凝土大坝，可以减小观测噪音和随机激励带来的干扰，鲁棒性好，并且对频率、振型的辨识精度较好，此外本方法在识别模态参数的计算效率高，可以更高效的为水利工程结构的抗震分析、健康诊断和震后结构损伤评估提供基础。
以下将结合附图及实施例对本发明做进一步详细说明。
附图说明
图1为实施例1中的四自由度弹簧-阻尼器-质量系统。
图2为实施例1中的四自由度时变结构模拟工况1采用不同时变模态参数辨识方法前两阶自振频率识别结果的对比。
图3为实施例1中的四自由度时变结构模拟工况2采用不同时变模态参数辨识方法前两阶自振频率识别结果的对比。
图4为实施例1中的四自由度时变结构模拟工况i在不同信噪比情况下的前两阶自振频率识别结果。
图5为实施例1中的四自由度时变结构模拟工况ii在不同信噪比情况下的前两阶自振频率识别结果。
图6为实施例1中模拟工况ii下snr＝40db时模态振型的识别值和理论值之间的模态置信因子(mac)值。
图7为pacoima拱坝在newhall强震中1、2、5、6、7和8径向观测通道的观测信息。
图8为采用newhall地震强震观测数据前n1＝400个样本绘制的奇异值谱图。
图9为用newhall地震数据得到的pacoima拱坝时变模态振型图。
图10为使用sanfernando地震数据得到的pacoima拱坝时变模态识别结果。
图11为使用chinohill地震数据得到的pacoima拱坝时变模态识别结果。
图12为使用newhall地震数据得到的pacoima拱坝时变模态识别结果。
图13为本发明的时变结构模态参数辨识方法的基本步骤。
具体实施方式
为进一步阐述本发明达成预定目的所采取的技术手段及功效，以下结合附图及实施例对本发明的具体实施方式详细说明如下。
实施例1：随机激励下的四自由度时变结构的模态参数辨识
使用如图1所示的4自由度的弹簧-质量-阻尼系统作为算例来证明本发明提出的时变模态参数辨识算法的识别精度、鲁棒性及计算效率。
其中系统的参数：m1＝m2＝m3＝1.0kg,m4＝0.9kg,k3＝7000n/m,k2＝k4＝8000n/m,c1＝c2＝0.6n·s2/m,c3＝c4＝0.55n·s2/m，时变的刚度系数k1变化情况如下所示，
模拟工况i：k1为突变性变化，
模拟工况ii：k1为线性连续变化，
k1＝6000+125t,(单位:n/m)。
系统的初始状态设为零，基础环境激励ag(t)采用高斯白噪音激励来模拟。结构振动响应的计算采用runge-kutta算法来实现。
本实施例中，基于广义子空间溯踪的时变结构模态参数辨识方法的具体实现步骤如下：
步骤(1)，根据香农采样定理，确定模拟加速度观测的采样频率为1000hz，振动信号的观测时间为16s，采集动力学响应信号得到离散的观测数据组成的向量y0(k)。
步骤(2)，根据本算例的工况在0.95到1.0之间调整遗忘因子β的取值，设定遗忘因子β＝0.995时模态参数识别精度较高，因此设定β＝0.995；初始数据长度n1取值越大，参数识别精度有所提高，同时计算时间也会延长，因此，本工况中，综合考虑计算效率与精度之后选取，选取n1＝400；采用初始的n1＝400个样本，计算初始的h0(k)并进行奇异值分解，根据特征值分解的结果绘制的奇异值谱，奇异值谱图中第一个转折处对应的奇异值数为8，因此，结构的系统阶数可以选为n＝8；在满足p×2l≥n的前提下(其中，n为系统阶次，l为时延向量的维数)，本实施例中选择p＝50。
步骤(3)，针对步骤(1)所采集的加速度响应信号中400个样本进行分析，采用希尔伯特变换获得相位偏移90度的虚拟响应y90(k)，构造以下的时延向量
并建立由该时延向量组成的矩阵
式(2)中，p为时延，满足p×2l≥n；l表示时延向量的维数，上标t表示转置，n表示系统阶次，
令k＝n1，根据式(2)得到时延向量组成的矩阵y(k)，然后，计算其相关函数矩阵h0(k)，得到的hankel矩阵h0(k)是一个方阵：
对相关函数矩阵h0(k)进行以下的特征值分解；
式(4)中，矩阵q1(k)包含了h0(k)的前n个特征向量，矩阵q1(k)是振动信号的主子空间，代表不同阶模态的贡献，q2(k)是剩余特征值对应的特征向量，δ1(k)和δ2(k)分别是由前n个主特征值和剩余特征值组成的对角矩阵，n表示系统阶次，上标t表示转置；
对正交矩阵π(k)、加速模态度响应的相关函数矩阵
π(k)＝q1(k)；(5)
σn(k)2＝h0(k)的最小特征值(7)。
步骤(4)，令k＝k+1，在tk时刻，获得新的观测数据y0(k)后，通过引入遗忘因子β对hankel矩阵h0(k)进行在线更新，hankel矩阵h0(k)在t＝tk时刻可以更新为
h0(k)＝βh0(k-1)+y(k)y(k)t(8)。
步骤(5)，为了计算加速模态度响应
对于广义子空间溯踪算法，首先引入以下式(9)的近似：
式(9)中，上标h表示复共轭转置，
这时，增广的(n+1)维列空间span(v(k))的标准正交基π(k)可表示为：
π(k)＝[π(k-1),u(k)](10)
式(10)中，u(k)是y⊥(k)的单位范数矢量，
定义矩阵
式(11)中，
式(12)中，q是向量qn(k)的前n个元素组成的向量，rn(k)是向量qn(k)的最后一个元素；
这时可以按照式(13)～(22)计算以下的参数b(k)、s(k)、c(k)、λm(k)和σn(k)2，向量e′(k)和q′(k)，以及矩阵π(k)和c2(k)，以实现对正交矩阵π(k)的在线更新：
b(k)＝1/(1+|rn(k)|)(13)
s(k)＝rn(k)/(|rn(k)|σ(k))(14)
π(k)＝π(k-1)-e′(k)qh(16)
σn(k+1)2＝min{λm(k),σn(k)2}(22)。
步骤(6)，按照步骤(5)的过程对正交矩阵π(k)进行更新后，令矩阵ο＝π(k)，在第k个时间间隔中，离散状态空间矩阵ak和观测矩阵gk使用下面的公式(23)和(24)来进行估计：
gk＝o(1:2l,:)(23)
式(23)和(24)中ο由矩阵o除掉最后的2l行元素后剩余的元素组成，矩阵
得到离散状态空间矩阵ak之后，固有频率和阻尼比使用其特征值进行提取，时变结构的固有模态频率fi从下式获得：
式(25)中，n0是模态阶数，对复模态系统，系统阶数n＝2n0，λi(i＝1,2,…,n0)是ak的特征值，re{·}和im{·}分别表示λi的实部和虚部；
时变结构的阻尼比ξi(i＝1,2,…,n0)使用下式(26)获得：
这时l个测量自由度对应的振型的分量φl，可以表达为
φl＝gkφ(27)
式(27)中，矩阵φ＝[φ1,φ2,...,φn]由ak的n个特征向量组成，φn是ak的第n个特征向量。
步骤(7)，继续运行上述计算步骤(4)～步骤(6)，直至所有的振动观测数据都已分析完成。
上述模态参数辨识过程的流程图如图13所示。
在信噪比snr＝40db水平下，采用增强的辅助变量子空间投影递归随机子空间识别(eiv-past-rssi)和基于广义子空间溯踪递归随机子空间识别(gyast-rssi)方法对模拟工况i、ii进行模态识别，其自振频率识别结果见图2和图3。图中黑色实线为常规eiv-past-rssi识别方法的结果，灰色虚线为采用本发明的gyast-rssi识别方法得到的结果；对比可得，采用eiv-past-rssi和gyast-rssi两种方法识别得到的结果差异性较小，但gyast-rssi的识别结果更加稳定。
为了证明该时变模态识别方法的鲁棒性，在不同信噪比(snr)水平下，采用时变模态识别方法gyast-rssi模拟的振动响应信号进行识别。对于工况i和ii，识别到的前两阶自振频率结果分别见图4和图5。图4和图5中细实线代表自振频率的理论值，灰色粗线为识别结果。从图4和图5中可以看出，对于工况i和ii，随着观测噪音的增强(信噪比减小，即由图(b)到图(d))，自振频率的识别精度依然很高，与理论值基本吻合，证明本文提出的基于gyast-rssi的时变结构模态识别方法的鲁棒性很强。
对于模拟工况ii，根据模态振型的识别值和理论值计算的模态置信因子(mac)见图6。从图6中可以看，四阶模态的mac值均大于0.99，说明gyast-rssi法对振型的识别精度比较高，采用该方法可以有效地跟踪时变结构模态振型的时变特性。
实施例2：pacoima拱坝在不同地震中的模态参数识别
将本发明的模态参数识别方法应用于实际的水利工程中，已验证本方法的实用性。
pacoima拱坝是美国南加州洛杉矶的一座高113m的拱坝。由于pacoima拱坝离著名的sanfernando地震危险区比较近，经历过多次不同震级地震，其中圣费尔南多地震(1971年)和北岭大地震(1994年)为该坝经历的两次较为严重的强震，由于圣费尔南多地震，坝体左坝肩支墩的构造缝张开，对大坝造成损伤。震后对坝体进行了修复工作，并且于1977年在坝体及河谷位置共布置了9个加速度传感器，共17个观测通道，自从1994年北岭大地震后，pacoima拱坝更新了该数字化强震观测系统，先后记录了3次地震响应，如表1所示。
表1：
步骤(1)，根据香农采样定理设定地震观测数据的采样频率为200hz，截止频率设定为50hz，对地震加速度响应信号进行采集得到离散的观测数据组成的向量y0(k)。强震观测系统记录的pacoima拱坝的1、2、5、6、7和8径向通道newhall地震加速度响应如图7所示。
步骤(2)，根据本算例的工况在0.95到1.0之间调整遗忘因子的取值，设定遗忘因子β＝0.995时，模态参数识别精度较高，因此设定β＝0.995；初始数据长度n1取值越大，参数识别精度有所提高，同时计算时间也会延长，因此，本工况中，综合考虑计算效率与精度之后选取，选取n1＝400；采用初始的n1＝400个样本，计算初始的h0(k)并进行奇异值分解，根据特征值分解的结果绘制的奇异值谱如图8所示，从图8中可以看出，奇异值谱图中第一个转折处对应的奇异值数为10，因此，大坝结构的系统阶数可以选为n＝10，模态阶数n0＝5；在满足p×2l≥n的前提下，本实施例中选择p＝100。
步骤(3)，将步骤(1)采集的地震加速度响应信号，用希尔伯特变换获得相位偏移90度的虚拟响应y90(k)，构造以下的数据向量
并建立由该时延向量组成的矩阵
式(2)中，p为时延，满足p×2l≥n，l表示时延向量的维数，上标t表示转置，n表示系统阶次，
令k＝n1，根据式(2)得到时延向量组成的矩阵y(k)，然后，计算其相关函数矩阵h0(k)，得到的hankel矩阵h0(k)是一个方阵：
对相关函数矩阵h0(k)进行以下的特征值分解：
式(4)中，矩阵q1(k)包含了h0(k)的前n个特征向量，矩阵q1(k)是振动信号的主子空间，代表不同阶模态的贡献，q2(k)是剩余特征值对应的特征向量，δ1(k)和δ2(k)分别是由前n个主特征值和剩余特征值组成的对角矩阵，n表示系统阶次，上标t表示转置；
对正交矩阵π(k)、加速模态度响应的相关函数矩阵
π(k)＝q1(k)；(5)
σn(k)2＝h0(k)的最小特征值(7)。
步骤(4)，令k＝k+1，在tk时刻，获得新的观测数据y0(k)后，通过引入遗忘因子β对hankel矩阵h0(k)进行在线更新，hankel矩阵h0(k)在t＝tk时刻可以更新为
h0(k)＝βh0(k-1)+y(k)y(k)t(8)。
步骤(5)，为了计算加速模态度响应
对于广义子空间溯踪算法，首先引入以下式(9)的近似：
式(9)中，上标h表示复共轭转置，
这时，增广的(n+1)维列空间span(v(k))的标准正交基π(k)可表示为：
π(k)＝[π(k-1),u(k)](10)
式(10)中，u(k)是y⊥(k)的单位范数矢量，
定义矩阵
式(11)中，
式(12)中，q是向量qn(k)的前n个元素组成的向量，rn(k)是向量qn(k)的最后一个元素；
这时可以按照式(13)～(22)计算以下的参数b(k)、s(k)、c(k)、λm(k)和σn(k)2，向量e′(k)和q′(k)，以及矩阵π(k)、c2(k)，以实现对正交矩阵π(k)的在线更新：
b(k)＝1/(1+|rn(k)|)(13)
s(k)＝rn(k)/(|rn(k)|σ(k))(14)
π(k)＝π(k-1)-e′(k)qh(16)
σn(k+1)2＝min{λm(k),σn(k)2}(22)。
步骤(6)，按照上述步骤(5)的过程对正交矩阵π(k)进行更新后，令矩阵ο＝π(k)，在第k个时间间隔中，离散状态空间矩阵ak和观测矩阵gk使用下面的公式(23)和(24)来进行估计：
gk＝o(1:2l,:)(23)
式(23)和(24)中ο由矩阵o除掉最后的2l行元素后剩余的元素组成，矩阵
得到离散状态空间矩阵ak之后，固有频率和阻尼比使用其特征值进行提取，时变结构的固有模态频率fi从下式获得：
式(25)中，n0是模态阶数，对复模态系统，系统阶数n＝2n0，λi(i＝1,2,…,n0)是ak的特征值，re{·}和im{·}分别表示λi的实部和虚部；
时变结构的阻尼比ξi(i＝1,2,…,n0)使用下式(26)获得：
这时l个测量自由度对应的振型的分量φl，可以表达为
φl＝gkφ(27)
式(27)中，矩阵φ＝[φ1,φ2,...,φn]由ak的n个特征向量组成，φn是ak的第n个特征向量。
步骤(7)，继续运行上述计算步骤(4)～步骤(6)，直至所有的振动观测数据都已分析完成。
上述模态参数辨识过程的流程图如图13所示。
对于三个典型的时刻，t＝20s,25s,30s和35s，使用newhall强震观测数据而识别得到的模态振型图结果见图9。图9中虚线代表的是原拱圈，实线代表各阶模态振型。
从图9中可以看出，在地震过程中的大坝模态振型的时变特性。
使用时变模态识别方法gyast-rssi对结构的自振频率和阻尼进行识别后，得到的识别结果见图10～图12。从图10～12中可以看出，随着地震响应幅值的变化，识别出的结构自振频率也相应的发生变化，说明大坝-基础-水库系统在地震激励下确实显示出了一定的时变特性。在地震响应幅值小的时段，自振频率识别值相对稳定；在地震响应幅值变化较大的时段，自振频率识别值有所减小。这是由于地震幅值较大时段，拱坝和基础上的一些接缝或者断层可能会张开，减小了结构刚度。
因此，实施例2证明了本发明的时变结构模态参数辨识方法用于水利工程时变结构，如混凝土大坝时，可以减小观测噪音和随机激励带来的干扰，鲁棒性好，并且对频率、振型的辨识精度较好，计算效率高，可以更高效的为水利工程结构的抗震分析、健康诊断和震后结构损伤评估提供基础。
以上内容是结合具体的优选实施方式对本发明所作的进一步详细说明，不能认定本发明的具体实施只局限于这些说明。对于本发明所属技术领域的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干简单推演或替换，都应当视为属于本发明的保护范围。
- 还没有人留言评论。精彩留言会获得点赞！