吳 杰， 徐洪俊， 張其林

(同濟(jì)大學(xué) 土木工程學(xué)院，上海 200092)

不確定度是評(píng)價(jià)測(cè)量結(jié)果質(zhì)量高低的重要指標(biāo)。隨著科學(xué)技術(shù)和生產(chǎn)的發(fā)展，不確定度分析受到越來(lái)越高的重視[1]。但在土木工程中，大部分的模態(tài)參數(shù)識(shí)別方法沒(méi)有考慮識(shí)別結(jié)果的不確定度問(wèn)題。Yuen等[2-3]提出了貝葉斯快速傅里葉變換(Bayesian Fast Fourier Transform，Bayesian FFT)方法，通過(guò)對(duì)原始時(shí)域數(shù)據(jù)進(jìn)行變換，識(shí)別結(jié)構(gòu)的模態(tài)特性并求解參數(shù)的不確定性。但該方法識(shí)別過(guò)程存在以下三個(gè)問(wèn)題：①由于后驗(yàn)的概率密度函數(shù)(Probability Density Function, PDF)是模態(tài)參數(shù)的非線(xiàn)性函數(shù)，因此模態(tài)參數(shù)的最佳估計(jì)(Most Probable Value, MPV)需要求解一個(gè)多維的數(shù)值優(yōu)化問(wèn)題，計(jì)算過(guò)程復(fù)雜；②Hessian矩陣需要通過(guò)有限差分法來(lái)確定，導(dǎo)致計(jì)算效率低下；③由于模態(tài)參數(shù)的數(shù)量隨著測(cè)試自由度數(shù)量的增加而增加，目標(biāo)函數(shù)的求解涉及到一個(gè)病態(tài)矩陣的求逆問(wèn)題。這些問(wèn)題的存在嚴(yán)重制約了該方法的應(yīng)用。

Au[4]針對(duì)選定的某一共振頻率帶的單個(gè)模態(tài)，提出了快速貝葉斯FFT(Fast Bayesian FFT)方法，該方法針對(duì)某一共振頻率帶的單個(gè)模態(tài)，不需要考慮所計(jì)算的自由度數(shù)量，模態(tài)參數(shù)的最佳估計(jì)可通過(guò)一個(gè)四維的數(shù)值優(yōu)化問(wèn)題來(lái)確定，而且Hessian矩陣可以通過(guò)對(duì)數(shù)似然函數(shù)關(guān)于模態(tài)參數(shù)的二階導(dǎo)數(shù)快速準(zhǔn)確地計(jì)算。

Fast Bayesian FFT方法雖然解決了模態(tài)參數(shù)的MPV和Hessian病態(tài)矩陣的求解問(wèn)題，但并沒(méi)有進(jìn)一步對(duì)模態(tài)參數(shù)的后驗(yàn)概率分布進(jìn)行探討[5-9]。鑒于此，本文首先基于貝葉斯理論和結(jié)構(gòu)響應(yīng)概率模型，提出了一種將多測(cè)點(diǎn)多自由度信號(hào)轉(zhuǎn)換成單測(cè)點(diǎn)單自由度信號(hào)的Bayesian TDD (Time Domain Decomposition)方法。然后，對(duì)單測(cè)點(diǎn)單自由度響應(yīng)下的貝葉斯FFT(Bayesian FFT)方法進(jìn)行了推導(dǎo)，并利用蒙特卡洛抽樣分析得出各個(gè)參數(shù)的后驗(yàn)概率分布、模態(tài)參數(shù)MPV和不確定性評(píng)估，從而避免了Hessian矩陣的求解運(yùn)算。最后，提出Bayesian TDD-FFT方法，并利用該方法對(duì)上海中心大廈實(shí)測(cè)數(shù)據(jù)進(jìn)行分析，探討了評(píng)估模態(tài)參數(shù)的后驗(yàn)概率分布及聯(lián)合分布問(wèn)題，同時(shí)將識(shí)別參數(shù)結(jié)果與Fast Bayesian FFT方法進(jìn)行了對(duì)比分析。

1 理論推導(dǎo)

1.1 Bayesian TDD方法

對(duì)于線(xiàn)性結(jié)構(gòu)，一般假設(shè)輸入為平穩(wěn)的零均值的Gauss過(guò)程，則結(jié)構(gòu)的輸出為0均值平穩(wěn)Gauss過(guò)程，即：

Xn×1～N(0,C)

(1)

式中:C為Xn×1響應(yīng)輸出的協(xié)方差矩陣。

根據(jù)多維高斯分布性質(zhì)[10], 對(duì)于n維的高斯分布Xn×1～N(0,C)，必存在正交陣Φn×n，使得：

X～N(0,ΦΛΦT)

(2)

Yn×1=ΦTX～N(0,Λ),其中Λ=diag(λi)

(3)

即Y的各個(gè)變量獨(dú)立。令θ={C},P(θ)為病態(tài)的均勻分布：P(θ)=1，滿(mǎn)足

P(θ|X)∝P(X|θ)P(θ)

(4)

(5)

(6)

將P(X|θ)寫(xiě)成對(duì)數(shù)形式，則有：

P(X|θ)=exp[-f(θ)]

(7)

其中，

(8)

由于ΦTCΦ=Λ,Λ=diag(λi)，則：

(9)

(10)

要使得P(θ|X)取得MPV值，需f(λ1,…,λi,…,λn)取最小值，則：

(11)

(12)

(13)

在實(shí)際結(jié)構(gòu)中，往往只關(guān)心所測(cè)結(jié)構(gòu)的前k階振型，測(cè)點(diǎn)的維數(shù)n>k，由振型疊加原理可得結(jié)構(gòu)的振動(dòng)響應(yīng)為：

(14)

1.2 單測(cè)點(diǎn)單自由度的Bayesian FFT

對(duì)于單測(cè)點(diǎn)單自由度振動(dòng)響應(yīng)，Bayesian FFT理論如下：

Yk=Fk+iGk=

(15)

考慮到采樣通道的電壓偏置和傅里葉變換的共軛等原因，僅k=1,2,…,Nq(Nq=int[N/2]+1)所對(duì)應(yīng)的FFT數(shù)據(jù)被應(yīng)用于結(jié)構(gòu)的模態(tài)參數(shù)識(shí)別，其所對(duì)應(yīng)的頻率值為fk=(k-1)/NΔt。

有關(guān)模態(tài)參數(shù)θ的后驗(yàn)概率密度函數(shù)P(θ|{Zk})正比于似然函數(shù)P({Zk}|θ)：

P(θ|{Zk})∝P({Zk}|θ)=

(16)

式中：Zk=[FkGk]T∈R2，k=1,2,3,…,Nq；Ck為矩陣Zk的協(xié)方差矩陣，其表達(dá)式為：

(17)

式中：σ2為預(yù)測(cè)誤差的功率譜密度；I2為2階單位矩陣；Hk為模態(tài)響應(yīng)自功率譜，表達(dá)式為：

(18)

(19)

式中：βk=f/fk為頻率比，f和ζ分別是單階模態(tài)信號(hào)特征頻率和阻尼比；S是單階模態(tài)響應(yīng)的自功率譜。式(17)可寫(xiě)成：

Ck=[(SDk+σ2)/2]I2

(20)

detCk=2-2(SDk+σ2)2

(21)

若將后驗(yàn)的PDF進(jìn)一步寫(xiě)成關(guān)于對(duì)數(shù)似然函數(shù)L(θ)的表達(dá)式，則有：

P(θ|{Zk})∝exp[-L(θ)]

(22)

式中：

(23)

將Ck和Zk代入式(23)，則有：

L(θ)=-(Nq-1)ln 2+

(24)

(25)

f、ζ、S和σ2的MPV值可通過(guò)對(duì)式(25)的無(wú)約束優(yōu)化直接求得。求得的f、ζ、S和σ2的MPV值后，利用蒙特卡羅方法對(duì)式(22)進(jìn)行采樣分析，得到f、ζ、S和σ2的后驗(yàn)概率分布，進(jìn)行不確定度分析。

2 計(jì)算流程

上文提出的Bayesian TDD-FFT方法的計(jì)算流程見(jiàn)圖1所示。

圖1 Bayesian TDD-FFT方法計(jì)算流程圖

3 數(shù)值驗(yàn)證

3.1 基本參數(shù)

本節(jié)通過(guò)數(shù)值模擬對(duì)Bayesian TDD-FFT模態(tài)參數(shù)識(shí)別方法進(jìn)行驗(yàn)證。利用中心差分法計(jì)算懸臂梁在隨機(jī)荷載作下的結(jié)構(gòu)響應(yīng)數(shù)據(jù)。懸臂梁詳細(xì)參數(shù)如下：長(zhǎng)度12 m，熱軋工字鋼I10，面積為14.345 cm2，慣性矩為245.0 cm4，Q345鋼材，楊氏模量E=206 GPa，泊松比為0.3，密度為7 850 kg/m3，取前3階模態(tài)進(jìn)行分析，每階模態(tài)的阻尼比均取0.5%，前三階頻率分別為0.806 8 Hz，5.056 8 Hz，14.160 6 Hz，荷載由隨機(jī)白噪聲生成，并為結(jié)構(gòu)響應(yīng)添加隨機(jī)白噪聲。沿梁長(zhǎng)度方向均勻設(shè)置12個(gè)加速度測(cè)點(diǎn)，采集頻率50 Hz，時(shí)長(zhǎng)1 200 s，加載點(diǎn)位于懸臂端，如圖2所示。圖3顯示了測(cè)點(diǎn)Ac3的加速度時(shí)程曲線(xiàn)。

圖2 懸臂梁測(cè)點(diǎn)布置

3.2 Bayesian TDD-FFT參數(shù)識(shí)別

圖4顯示了測(cè)點(diǎn)Ac3的加速度功率譜。利用AMD(Analytical Mode Decomposition)[12]法對(duì)原始數(shù)據(jù)進(jìn)行濾波處理，帶寬分別取0～2 Hz、4～6 Hz、13～15 Hz。利用Bayesian TDD法計(jì)算得到懸臂梁的前3階模態(tài)振型如圖5所示，從圖中可看出，計(jì)算值與理論值吻合很好，表明Bayesian TDD方法能有效地識(shí)別結(jié)構(gòu)振型。

圖3 測(cè)點(diǎn)Ac3加速度時(shí)程

圖4 測(cè)點(diǎn)Ac3加速度功率譜

圖5 Bayesian TDD-FFT法識(shí)別出的前3階振型

Fig.5 The first three modes identified by Bayesian TDD-FFT method

表1列出了部分模態(tài)參數(shù)的識(shí)別結(jié)果，其中，頻率的識(shí)別值與理論值0.806 8 Hz、5.056 8 Hz、14.160 6 Hz很接近，阻尼比的識(shí)別值與理論值0.005也基本一致，表明Bayesian TDD-FFT方法能有效地識(shí)別結(jié)構(gòu)的動(dòng)力參數(shù)。從變異系數(shù)(Coefficient of Variation，COV)可以看出，阻尼比的變異系數(shù)是最高的，可見(jiàn)在模態(tài)識(shí)別中阻尼比的值不易準(zhǔn)確測(cè)定。

4 工程應(yīng)用

4.1 上海中心加速度響應(yīng)采集系統(tǒng)

為了采集自然激勵(lì)下的加速度數(shù)據(jù)，上海中心大廈共沿樓層高度方向共安裝了71個(gè)Lance加速度傳感器。加速度傳感器的詳細(xì)參數(shù),見(jiàn)表2。

表1 Bayesian TDD-FFT模態(tài)參數(shù)識(shí)別結(jié)果

Tab.1 Modal parameters identification based on Bayesian TDD-FFT method

階次參數(shù)頻率阻尼比激勵(lì)譜密度誤差譜密度第一階MVP0.807 40.004 63.823 71.798 5E-2標(biāo)準(zhǔn)差5.368 4E-34.103 2E-41.727 2E-31.490 5E-4COV/%0.664 98.920 04.517 1E-20.828 7第二階MVP5.050 20.004 92.611 31.262 0E-2標(biāo)準(zhǔn)差0.039 44.790 4E-46.072 6E-41.503 4E-4COV/%0.780 29.855 52.325 5E-21.191 3第三階MVP14.167 80.004 82.168 11.940 8E-3標(biāo)準(zhǔn)差0.064 763.843 8E-45.444 2E-43.655 5E-5COV/%0.457 18.007 12.511 0E-21.883 5注:變異系數(shù)COV=標(biāo)準(zhǔn)差/MPV

表2 加速度傳感器參數(shù)

4.2 良態(tài)風(fēng)作用下模態(tài)參數(shù)識(shí)別

2018年4月12日～14日上海浦東風(fēng)速達(dá)到4～5級(jí)，尤其是4月13日的瞬時(shí)最大風(fēng)速達(dá)到30 m/s，實(shí)時(shí)風(fēng)速見(jiàn)圖6?？紤]到識(shí)別效果和儀器的采樣精度，取風(fēng)速較大時(shí)間段(13日7時(shí)～8時(shí))的加速度數(shù)據(jù)進(jìn)行模態(tài)參數(shù)識(shí)別?？紤]到篇幅，本文只選取21層、36層、51層、67層、83層、100層、117層、124層的X(東西)向加速度進(jìn)行識(shí)別。圖7是各層加速度的樣本(600 s)和117層加速度功率譜密度。

圖6 2018年4月13日實(shí)時(shí)風(fēng)速

利用AMD(Analytical Mode Decomposition)法對(duì)原始數(shù)據(jù)進(jìn)行濾波處理，帶寬分別取[0，0.17]、[0.17，0.25]。對(duì)濾波后的數(shù)據(jù)應(yīng)用Bayesian TDD-FFT方法進(jìn)行處理，計(jì)算結(jié)果見(jiàn)圖8～14和表3所示。其中，圖8和圖9顯示了前兩階模態(tài)和單自由度響應(yīng)概率分布。前兩階模態(tài)與Fast Bayesian FFT求解結(jié)果相似。單自由度響應(yīng)的概率密度分布與Gauss分布十分吻合，驗(yàn)證了Bayesian TDD假定和推論的正確性。

表3 Bayesian TDD-FFT前兩階模態(tài)參數(shù)

Tab.3 The first two order modal parameters based on Bayesian TDD-FFT method

(a) 不同樓層600 s加速度樣本

(b) 117層加速度功率譜

Fig.7 Acceleration sample during the 600 s and power spectral density of 117th floor acceleration

X方向前兩階模態(tài)參數(shù)抽樣樣本如圖10所示，MPV值見(jiàn)表3。由于抽樣起始點(diǎn)是MPV值，因此抽樣樣本收斂在MPV值附近。從圖11和圖12可以看出，頻率和阻尼比的后驗(yàn)概率密度分布與Gauss分布相吻合，而激勵(lì)譜密度和誤差譜密度的后驗(yàn)概率密度分布離散且與Gauss分布相差較大。

圖13和圖14為前兩階的模態(tài)參數(shù)后驗(yàn)聯(lián)合分布直方圖，從圖中可看出，頻率與阻尼比沒(méi)有相關(guān)性，而激勵(lì)譜密度和預(yù)測(cè)誤差譜密度呈近似的線(xiàn)性相關(guān)性：預(yù)測(cè)誤差譜密度與激勵(lì)譜密度呈反比例關(guān)系，即噪聲譜密度隨著激勵(lì)譜密度的增加而降低，這與信號(hào)采集時(shí)的常規(guī)認(rèn)知相一致。

(a) X方向第一階陣型

(b) X方向第一階模態(tài)響應(yīng)概率分布

Fig.8 The 1stmode shape and the probability distribution of the 1stmode response inXdirection

(a) X方向第二階陣型

(b) X方向第二階模態(tài)響應(yīng)概率分布

Fig.9 The 2ndmode shape and the probability distribution of the 2ndmode response inXdirection

(a) 第一階模態(tài)參數(shù)抽樣樣本(b) 第二階模態(tài)參數(shù)抽樣樣本

圖10 X向第一階模態(tài)和第二階模態(tài)參數(shù)抽樣樣本

圖11 第一階模態(tài)參數(shù)后驗(yàn)概率密度分布

Fig.11 The posterior probability density distribution of the 1stmodal parameters

(a) 頻率后驗(yàn)概率分布圖(b) 阻尼比的后驗(yàn)概率分布圖(c) 激勵(lì)譜的后驗(yàn)概率分布圖(d) 誤差譜后驗(yàn)概率分布圖

圖12 第二階模態(tài)參數(shù)后驗(yàn)概率密度分布

Fig.12 The posterior probability density distribution of the 2ndmodal parameters

表4為應(yīng)用Fast Bayesian FFT方法計(jì)算得到的前兩階模態(tài)結(jié)果，F(xiàn)ast Bayesian FFT法采用的濾波帶寬與Bayesian TDD-FFT法相同。從表3和表4可看出：Bayesian TDD-FFT法所求得的頻率和阻尼與Fast Bayesian FFT法接近，但前者的變異系數(shù)更低，結(jié)果更穩(wěn)定。Bayesian TDD-FFT法所求激勵(lì)譜密度大于Fast Bayesian FFT法所求值，而誤差譜密度小于后者，這是由于在Bayesian TDD法處理時(shí)舍棄了部分噪聲分量。表中變異系數(shù)表明，兩種方法求解出的頻率、激勵(lì)譜密度和誤差譜密度的后驗(yàn)不確定性小，而阻尼比的不確定性大。

表4 Fast Bayesian FFT前兩階模態(tài)參數(shù)

Tab.4 The first two order modal parameters based on Bayesian FFT method

階次參數(shù)頻率阻尼比激勵(lì)譜密度誤差譜密度第一階MVP0.1060.0141.193E-066.114E-05標(biāo)準(zhǔn)差6.546E-42.206E-39.278E-81.190E-8COV/%0.6115.927.771.947 2E-02第二階MVP0.2030.0131.226E-071.241E-06標(biāo)準(zhǔn)差3.854E-41.446E-39.526E-95.479E-12COV/%0.1910.987.764.414 9E-04注:變異系數(shù)COV=標(biāo)準(zhǔn)差/MPV

表5列出了SSI(Stochastic Subspace Identification，隨機(jī)子空間)法和Bayesian TDD-FFT法識(shí)別的模態(tài)頻率和阻尼比的對(duì)比結(jié)果，從表中可看出，兩種方法的頻率識(shí)別結(jié)果基本相同，但阻尼比的識(shí)別結(jié)果存在一定差異，進(jìn)一步說(shuō)明了阻尼比識(shí)別的不確定性相對(duì)較大。相比SSI方法，Bayesian TDD-FFT法的優(yōu)勢(shì)在于能夠獲得參數(shù)的后驗(yàn)概率分布，同時(shí)可以對(duì)參數(shù)的不確定性進(jìn)行評(píng)估。

表5 Bayesian TDD-FFT和SSI方法參數(shù)識(shí)別結(jié)果對(duì)比

Tab.5 Result comparison of parameter identification using SSI and Bayesian TDD-FFT method

階次Bayesian TDD-FFTSSI頻率/Hz阻尼比頻率/Hz阻尼比第一階 0.1060.012 10.1060.014 2第二階 0.2030.015 60.2040.017 6

5 結(jié) 論

Bayesian TDD法將多測(cè)點(diǎn)多自由度信號(hào)轉(zhuǎn)換成單測(cè)點(diǎn)單自由度信號(hào)，有效地降低了監(jiān)測(cè)信號(hào)處理的復(fù)雜程度。Bayesian TDD和蒙特卡羅方法的聯(lián)合應(yīng)用,避免了Bayesian FFT法中病態(tài)Hessian矩陣的求逆運(yùn)算，能夠有效地對(duì)結(jié)構(gòu)模態(tài)參數(shù)進(jìn)行評(píng)估。實(shí)際工程應(yīng)用結(jié)果表明：

(a) 頻率與阻尼比后驗(yàn)聯(lián)合分布直方圖

(b) 頻率與阻尼比聯(lián)合分布熱力圖

(c) 激勵(lì)譜密度與誤差譜密度后驗(yàn)聯(lián)合分布直方圖

(d) 激勵(lì)譜密度與誤差譜密度聯(lián)合分布熱力圖

Fig.13 Posterior joint distribution histogram and hot map of the 1stmodal parameters

(1) Bayesian TDD-FFT法能夠有效地識(shí)別結(jié)構(gòu)的模態(tài)參數(shù)，利用蒙特卡羅抽樣分析可以獲得及驗(yàn)證模態(tài)參數(shù)的后驗(yàn)概率分布，此為Bayesian TDD-FFT法優(yōu)于Fast Bayesian FFT法之處；

(2) 對(duì)比分析表明：Bayesian TDD-FFT法雖然是一種數(shù)值解法，但仍能夠達(dá)到Fast Bayesian FFT法的解析解精度；相比SSI法，Bayesian TDD-FFT法的優(yōu)勢(shì)在于能夠獲得參數(shù)的后驗(yàn)概率分布，同時(shí)可以對(duì)參數(shù)的不確定性進(jìn)行評(píng)估；

(a) 頻率與阻尼比后驗(yàn)聯(lián)合分布直方圖

(b) 頻率與阻尼比聯(lián)合分布熱力圖

(c) 激勵(lì)譜密度與誤差譜密度后驗(yàn)聯(lián)合分布直方圖

(d) 激勵(lì)譜密度與誤差譜密度聯(lián)合分布熱力圖

Fig.14 Posterior joint distribution histogram and hot map of the 2ndmodal parameters

(3) 上海中心大廈的實(shí)測(cè)數(shù)據(jù)分析結(jié)果表明：頻率、激勵(lì)譜密度和誤差譜密度的后驗(yàn)不確定性小，而阻尼比的不確定性大；頻率和阻尼沒(méi)有相關(guān)性，而激勵(lì)譜密度和誤差譜密度存在近似的線(xiàn)性相關(guān)性，即反比例關(guān)系。