磁流變阻尼器MNS模型參數(shù)不確定性分析

彭長樂，陳 城，侯和濤

(山東大學(xué)土建與水利學(xué)院，山東，濟南 250000)

結(jié)構(gòu)控制的概念在土木工程領(lǐng)域提出以來，經(jīng)過多年來的研究發(fā)展，在減小結(jié)構(gòu)的動力荷載響應(yīng)方面顯示出了巨大的潛力[1-2]，其中磁流變阻尼器因為其低能耗、出力大、響應(yīng)快、性能可靠等優(yōu)點備受研究者們關(guān)注[3-4]。特別是自Lord公司研制出最大出力200 kN可用于土木工程領(lǐng)域的足尺阻尼器以來，磁流變阻尼器更成為結(jié)構(gòu)控制領(lǐng)域的熱點。要使磁流變阻尼器的控制性能得到充分的發(fā)揮，獲取一個準(zhǔn)確的數(shù)值模型是十分必要的。迄今為止，許多磁流變阻尼器的模型相繼被提出來，其中包括Bouc-Wen 模型[5]、Viscous-Dahl模型[6]，Hyperbolic Tangent模型[7]和Maxwell Nonlinear Slider (MNS)模型[8]等。這些模型的參數(shù)多以進化算法設(shè)計目標(biāo)函數(shù)，基于一系列的正弦位移識別試驗，通過優(yōu)化計算出最優(yōu)參數(shù)[9]。

然而由于在建模情況中所做的簡化假設(shè)，以及有限的觀測數(shù)據(jù)，傳統(tǒng)的確定性參數(shù)現(xiàn)象學(xué)模型很難考慮由于環(huán)境因素、模型誤差以及測量誤差帶來的不確定性，從而導(dǎo)致數(shù)值模型的出力結(jié)果與實驗結(jié)果存在不同程度的差異。這些差異可能使得數(shù)值模擬中的結(jié)構(gòu)在地震響應(yīng)下響應(yīng)結(jié)果與實際實驗不符合，也會一定程度影響磁流變阻尼器的控制算法效果。Caicedo等[10]證明使用模型參數(shù)的概率表征方式有助于更真實的模擬預(yù)測。

MNS模型利用Hershel-Bulkley粘塑性來描述阻尼器中發(fā)生的MR流體的后屈服非牛頓流體行為，即剪切稀化和增稠行為[8]。MNS模型能比較準(zhǔn)確描述阻尼器的位移與出力、速度與出力的非線性滯回行為。但由于MNS模型本身屬于現(xiàn)象學(xué)模型的一種，阻尼器的復(fù)雜非線性行為使得建模仍存在偏差，因此單一的確定性模型并不能完全滿足預(yù)測的要求。而不確定性參數(shù)能幫助克服參數(shù)優(yōu)化時出現(xiàn)“異參同效”的影響，最為有效地利用阻尼器模型的預(yù)測能力。

本文運用馬爾可夫蒙特卡洛方法對磁流變阻尼器的MNS模型進行了不確定性分析。通過選擇合適的似然函數(shù)生成模型參數(shù)的后驗分布，并對參數(shù)分布形式，參數(shù)之間的內(nèi)在相關(guān)性進行分析。在此基礎(chǔ)上進一步與確定性參數(shù)在正弦位移曲線實驗和實時混合模擬實驗下的出力誤差、耗能預(yù)測方面進行對比總結(jié)，分析不確定性模型的特點。

1 MNS確定性模型介紹

磁流變阻尼器的MNS模型由屈服前和屈服后兩部分組成[9], 如圖1所示。

1.1 屈服前模式

屈服前的阻尼器行為由Maxwell單元描述，由一個系數(shù)為c的阻尼元件與一個剛度為k的彈簧元件串聯(lián)而成。其阻尼力f為：

圖1 MNS模型示意圖Fig.1 Schematic diagram of MNS model

式中：y表示彈簧元件與阻尼元件的總位移；z表示阻尼元件的位移。當(dāng)阻尼器處于屈服前模式，阻尼元件與彈簧元件的總速度與實際阻尼器速度相等。在阻尼器初始狀態(tài)時y與x相等，直到阻尼器出現(xiàn)第1次屈服。在屈服后模式下，y與z需要通過式(1)的變式，并利用式(2)計算的屈服后出力f值持續(xù)更新。

當(dāng)最終模式從屈服后返回屈服前，這時的y=x+σ，并且在此后的屈服前計算中σ恒定不變。

1.2 屈服后模式

當(dāng)阻尼器行為模式從預(yù)屈服變?yōu)榍?，阻尼力的斜率發(fā)生改變，阻尼力通常隨著速度的增加而增加。在屈服后模式下，MNS模型通過運用非線性滑塊(nonlinear slider)來描述阻尼器出力與速度之間的關(guān)系，將正負向出力分開兩組表示，從而有利于提高阻尼器不對稱下的預(yù)測。正向阻尼器出力的數(shù)學(xué)表達式如下：

第1部分基于Herschel-Bulkley粘塑性理論，包括、a、b、n這4個參數(shù)。第2部分是與Herschel-Bulkley粘塑性理論定義的切線。負向阻尼器出力可以類似公式表達，因此同樣也包含4個負向出力參數(shù)。從加載、卸載路徑上來看，磁流變阻尼器的出力也存在差異，可以利用增加式(6)的方式模擬不同階段下的出力。

1.3 模式改變判別

當(dāng)阻尼器出力從Maxwell單元達到屈服后曲線，模型從屈服前轉(zhuǎn)換屈服后。數(shù)學(xué)上這種狀態(tài)可以表達為：

當(dāng)模型從屈服后回到屈服前，只需要速度滿足:式(8)就會出現(xiàn)：

式中，y可以利用式(1)的變形計算：

2 試驗數(shù)據(jù)

本文使用的試驗數(shù)據(jù)基于科羅拉多大學(xué)智能結(jié)構(gòu)技術(shù)實驗室完成的足尺磁流變阻尼器性能測試試驗[10]。由于篇幅的限制，本文主要考慮阻尼器在最大電流2.5A下，即為Passive-On模式下的不確定性研究。表1列出了本文選取的6組正弦位移加載實驗的頻率幅值，每組都在2.5 A的電流輸入下進行10次往復(fù)加載。

表1 正弦試驗參數(shù)Table 1 Sinusoidal test parameters

通過利用粒子群算法[9]，設(shè)計最小均方根誤差為目標(biāo)函數(shù)，經(jīng)過多次迭代，計算出模型與試驗相對誤差最小的確定性MNS模型參數(shù)值列于表2中。

表2 MNS模型參數(shù)信息統(tǒng)計結(jié)果Table 2 Statistics of MNS model parameter

3 馬爾可夫鏈蒙特卡洛方法

馬爾科夫鏈蒙特卡洛算法可以獲得參數(shù)后驗分布的一系列樣本，適用于非標(biāo)準(zhǔn)分布不獨立的多變量分布模擬[11]。其基本步驟如下：

1) 建立先驗分布：由于對于本文模型沒有任何先驗信息。因此本文的所有參數(shù)的先驗分布Prior都假設(shè)為在區(qū)間內(nèi)的均勻分布。

2) 初始點：合適的參數(shù)初始值可以減少達到收斂所需的計算時間，提高運算效率。本文中選擇已有文獻中的確定性參數(shù)作為初始點。

3) 計算似然值：每組參數(shù)值會生成對應(yīng)的阻尼力模型輸出，通過計算與試驗輸出的相對誤差，得到模型擬合程度的函數(shù)如下：

假設(shè)似然函數(shù)符合均值為0，標(biāo)準(zhǔn)差為ε的高斯分布。ε的選值直接影響后驗參數(shù)的范圍，過小的ε會使參數(shù)收斂到一個狹窄的區(qū)域，不確定性太小。過大的ε則會使得參數(shù)范圍太大，導(dǎo)致不確定分析的結(jié)果過度偏離試驗結(jié)果。

4) 計算后驗概率：

5) 生成新樣本：加上隨機增量會新生成一組參數(shù)，新參數(shù)計算的后驗如果滿足式(13)則認為可以接受：

其中，γ服從[0，1]的均勻分布。重復(fù)步驟3)~步驟5)，生成長度為N的馬爾科夫鏈。選擇收斂后的鏈長M為參數(shù)樣本進行不確定性分析。

4 參數(shù)分析結(jié)果

通過選取ε值為0.5%，設(shè)定抽樣長度為十萬，得到了相應(yīng)的后驗分布參數(shù)樣本。本文利用Gelman-Rubin原則[12]，通過分析多個馬爾可夫鏈之間的差異來評估MCMC收斂。通過比較每個模型參數(shù)的估計的鏈間和鏈內(nèi)方差，確定MCMC的結(jié)果在50000次后基本達到收斂。最終選擇收斂后的馬爾科夫鏈作為合理近似參數(shù)樣本。

圖2做出不確定參數(shù)分布與粒子群算法得到最優(yōu)確定性參數(shù)的對比。可以發(fā)現(xiàn)所有參數(shù)都為單峰分布。其中c與k呈現(xiàn)正態(tài)分布，其他參數(shù)則表現(xiàn)出對數(shù)正態(tài)分布的趨勢。從圖2中還可以發(fā)現(xiàn)MNS模型的確定性參數(shù)值都處在不確定性參數(shù)的分布之中，并且大部分都位于概率密度函數(shù)的峰值附近，卻不一定恰好出現(xiàn)在頂峰位置。這主要是因為參數(shù)之間出現(xiàn)一些相關(guān)性，導(dǎo)致最優(yōu)的參數(shù)組合不一定都處于各自參數(shù)的最可能出現(xiàn)的點，這也從側(cè)面反映了參數(shù)之間的相互影響會導(dǎo)致確定性參數(shù)的取值，單組參數(shù)并不能完全發(fā)揮出模型的預(yù)測能力。表2列出了對各參數(shù)的統(tǒng)計信息進行分析得到的前3階矩。

圖2 不確定性參數(shù)分布與確定性參數(shù)對比Fig.2 Comparison of probabilistic parameter distribution and deterministic parameters

圖3 MNS模型參數(shù)變異系數(shù)Fig.3 Coefficient of variation for MNS model parameters

圖3使用變異系數(shù)來具體定量分析參數(shù)本身的不確定性大小。從圖3中可以觀測到，在由馬爾可夫鏈蒙特卡洛方法得到的不確定性參數(shù)中，不同參數(shù)分布的變異系數(shù)不盡相同。對于有明顯物理意義的參數(shù)，描述阻尼器在未屈服狀態(tài)下的阻尼系數(shù)c與剛度系數(shù)k的變異系數(shù)明顯要小于其他參數(shù)。而用于擬合阻尼器在屈服后Hershel-Bulkley曲線的模型參數(shù)，則顯示出較強的不確定性。這說明對于有較強客觀物理意義的參數(shù)，其識別能力更強，不確定性更小。對于現(xiàn)象學(xué)模型中的參數(shù)，由于其本身沒有明顯物理意義，一般為擬合曲線而設(shè)計，所存在的參數(shù)組合可能性更多，不確定性更大。

為考慮參數(shù)相關(guān)性的影響，本文對參數(shù)組合間的相關(guān)系數(shù)矩陣進行分析計算，參數(shù)之間的某些較強相關(guān)性說明參數(shù)相互間是多余的，可以為模型簡化提供指導(dǎo)。圖4可視化出參數(shù)之間的相關(guān)性大小，從中明顯可以看出：描述屈服前的參數(shù)阻尼系數(shù)c與剛度系數(shù)k存在較強的相關(guān)性，為0.5左右；描述屈服后的模式a，b，n之間同樣存在很強相關(guān)性，大致為0.6~0.9，并且僅僅出現(xiàn)在相同正負的方向中；其余的參數(shù)之間相關(guān)性不強，幾乎在-0.2到0.2之間，可以基本忽略相關(guān)性而認為這些參數(shù)是相互獨立的。

圖4 MNS模型參數(shù)相關(guān)系數(shù)Fig.4 Correlation between MNS model parameters

5 基于不確定性模型的試驗預(yù)測

模型參數(shù)不確定性分析的最終目的是利用模型參數(shù)對真實物理模型的響應(yīng)進行有效的預(yù)測。本節(jié)分別從試驗的阻尼力輸出和耗能預(yù)測兩個指標(biāo)，對不確定性模型和確定性模型進行對比評估。以表1中25.4 mm，1 Hz正弦位移加載試驗為例，圖5給出了磁流變阻尼器阻尼力輸出時程比較，其中實線為不確定性參數(shù)模擬的95%預(yù)測區(qū)間以及預(yù)測均值，虛線為確定性模型預(yù)測曲線，點劃線為試驗值?？梢钥吹酱_定性參數(shù)預(yù)測曲線幾乎與不確定性曲線的均值預(yù)測重合。另外試驗觀測值大部分都處于95%的預(yù)測區(qū)間之內(nèi)，表明不確定性模型具有良好的預(yù)測能力。值得指出，試驗值仍有處于預(yù)測之外的值，這主要受限制于模型本身的精確程度，參數(shù)的不確定性分析，并不能完全補償模型誤差的影響。

圖5 磁流變阻尼器輸出力比較Fig.5 Comparison of force outputs

為定量分析，本文使用歸一化均方根(Root Mean Square, RMS)誤差作為評價模型參數(shù)優(yōu)劣的目標(biāo)函數(shù)。

式中：ERMS越小，模型預(yù)測與實際輸入的差距越??；當(dāng)ERMS=0時，表示模型與試驗結(jié)果完全一致。對于表1中的六組模型識別試驗，圖6比較了MNS 模型對應(yīng)確定性和不確定性參數(shù)的ERMS。

圖6 模型誤差分布比較Fig.6 Comparison of RMS error

從圖6可以發(fā)現(xiàn)，本文中得到的不確定性參數(shù)的誤差基本處在3%~4%之間，而確定性算法得到一組參數(shù)值的誤差接近不確定性模型的誤差最小值處。在參數(shù)識別過程中，不確定性模型適當(dāng)放寬了接受參數(shù)的范圍，而非僅僅接受最優(yōu)的一組參數(shù)。使用這些參數(shù)的組合可以在合理范圍內(nèi)拓寬模型的預(yù)測區(qū)間，從而提高模型的預(yù)測性能。

為了進一步檢驗對MNS模型參數(shù)不確定性分析的預(yù)測效果，本文采用了Chen等[13]進行的實時混合模擬試驗結(jié)果。該系列混合模擬試驗選用了1999年Chi-Chi地震在TCU105站點觀測到的地震運動加速度并縮放到Design Basis Earthquake (DBE)和Maximum Considered Earthquake (MCE)兩種情況。磁流變阻尼器在DBE和MCE下實時混合模擬試驗(RTHS)的最大位移分別在50 mm和70 mm左右，對應(yīng)的磁流變阻尼器速度-出力曲線如圖7和圖8所示。

對應(yīng)不同情況下阻尼器的出力情況，計算出的均方根誤差ERMS分布如圖9和圖10所示。對于這兩組混合模擬試驗結(jié)果，確定性參數(shù)的阻尼力輸出誤差不再處于不確定性參數(shù)模型的誤差最小點。這說明對于不同的實驗，最優(yōu)參數(shù)并不是唯一的，這也同時說明了不確定性分析的重要性。完全靠一組參數(shù)值是很難做到在實驗情況出現(xiàn)變化時仍保證較高預(yù)測精度的。而在控制不確定性大小合適的情況下，可以彌補確定性參數(shù)的不足，從而提高MNS模型對不同試驗的預(yù)測能力。

圖7 實時混合模擬試驗DBE下速度與出力的關(guān)系Fig.7 Velocity-force curve of RTHS with DBE

圖8 實時混合模擬試驗MCE下速度與出力的關(guān)系Fig.8 Velocity-force curve of RTHS with MCE

圖9 實時混合模擬試驗中震下ERMSFig.9 ERMS for real time hybrid simulation with DBE

圖10 實時混合模擬試驗大震下ERMSFig.10 ERMS of real time hybrid simulation test with MCE

從圖9和圖10進一步可以看出不確定模型在正弦位移加載方式下的誤差基本處于3%~4%之間，而在實時混合模擬試驗的位移響應(yīng)下的誤差增加至7%~11%。這說明了MNS模型在有規(guī)律正弦位移加載方式下的預(yù)測能力要強于模型對于實時混合模擬試驗的預(yù)測。這提醒我們在參數(shù)優(yōu)化過程中，過分依賴使用有規(guī)律的位移加載方法會低估模型在真實結(jié)構(gòu)位移中的誤差。

磁流變阻尼器作為建筑結(jié)構(gòu)中的控制裝置，主要起耗散能量的作用。因此在用模型預(yù)測與實際出力的誤差來衡量參數(shù)優(yōu)劣外，本文進一步用耗能預(yù)測來檢驗?zāi)Ｐ筒淮_定性分析的預(yù)測效果。

首先對正弦位移識別實驗本身進行耗能預(yù)測，圖11對不確定性參數(shù)預(yù)測的耗能分布與正弦位移

圖11 不同頻率、幅值下耗能預(yù)測(直方圖為不確定性模型耗能預(yù)測，虛線為確定性模型預(yù)測耗能，點劃線為試驗耗能)Fig.11 Energy dissipation prediction under different frequencies and amplitudes (The histogram is the uncertainty prediction, the dashed line is the deterministic prediction, and the dotted line is the experimental energy consumption)

實驗實際耗能預(yù)測結(jié)果進行對比。從圖中可以觀測到：使用不確定性參數(shù)模型得到的耗能分布基本上都處于近似正態(tài)分布；確定性模型預(yù)測一般都接近于不確定性模型預(yù)測的均值處。并且試驗得到的實際耗能都處于不確定性模型預(yù)測區(qū)間內(nèi)部，這說明使用這組不確定性參數(shù)組合預(yù)測這6組試驗是完全可行的。

圖12對圖7和圖8中實時混合模擬試驗進行耗能預(yù)測分析。確定性參數(shù)預(yù)測仍接近不確定性模型預(yù)測均值附近，并且不確定性參數(shù)的預(yù)測耗能區(qū)間依然將試驗耗能包絡(luò)其中。但是值得注意的是，相對于正弦位移識別試驗。實時混合模擬下的試驗耗能，處于耗能預(yù)測的更邊緣位置，這說明模型的耗能預(yù)測能力下降，但是仍處于5%~95%預(yù)測范圍區(qū)間之內(nèi)。這與使用ERMS分析得到的結(jié)論類似，即正弦位移下MSN模型預(yù)測比實時混合模擬下的結(jié)構(gòu)響應(yīng)出力要更準(zhǔn)確。

圖12 不同地震水平下磁流變阻尼器耗能預(yù)測(直方圖為不確定性模型耗能預(yù)測，虛線為確定性模型預(yù)測耗能，點劃線為試驗耗能)Fig.12 Energy dissipation prediction of MR damper for real-time hybrid simulation under difterent seismic levers (The histogram is the uncertainty prediction, the dashed line is the deterministic prediction, and the dotted line is the experimental energy consumption)

6 結(jié)論

本文對200 kN足尺磁流變阻尼器的MNS模型在2.5 A下的參數(shù)不確定性進行了分析。通過選擇合適的不確定性大小與似然函數(shù)，基于馬爾科夫鏈蒙特卡洛方法，生成了一系列參數(shù)組合值。并使用后驗參數(shù)分別對正弦位移試驗和實時混合模擬下的阻尼器出力進行了預(yù)測評估。根據(jù)分析結(jié)果得出以下結(jié)論。

(1) 從生成參數(shù)的分布與相關(guān)性中分析發(fā)現(xiàn)：最優(yōu)的一組參數(shù)值被包含在不確定性參數(shù)的分布中，并且大部分確定性參數(shù)都趨近于均值處；某些參數(shù)間存在強烈的相關(guān)性，這為去除多余參數(shù)，簡化模型提供了思路；參數(shù)的變異系數(shù)不盡相同，有實際物理意義的較小，更容易識別，而現(xiàn)象模型參數(shù)則存在更多可能性。

(2) 利用不確定性分析的參數(shù)分布預(yù)測試驗出力可以發(fā)現(xiàn)，大部分試驗出力處于不確定性模型的5%~95%預(yù)測區(qū)間以內(nèi)，但是不足以完全包絡(luò)。這意味著不確定性模型可以提高預(yù)測能力，但是并不可能完全彌補模型本身的誤差。

(3) 在對歸一化均方根誤差的分析中，選用的正弦位移識別試驗時，確定性值的誤差最小，而在混合模擬驗證試驗中，這組確定性值的出力預(yù)測卻并不處于誤差最小值點，并且在實時混合模擬驗證試驗中誤差明顯變大了。說明單一最優(yōu)參數(shù)的適應(yīng)性較弱，而使用不確定性參數(shù)則可以減小由于試驗情況改變引起的誤差增大。

(4) 最后通過對正弦識別試驗和實時混合模擬驗證試驗進行耗能預(yù)測。所有選取的試驗中，確定性模型預(yù)測處于不確定性模型預(yù)測區(qū)間的均值處，而阻尼器實際耗能都處于預(yù)測區(qū)間內(nèi)部，這說明不確定性模型對于阻尼器耗能的預(yù)測是可行的，可以為進一步預(yù)測建筑結(jié)構(gòu)在使用磁流變阻尼器時的位移響應(yīng)預(yù)測提供參考。