不確定性參數(shù)識(shí)別的區(qū)間響應(yīng)面模型修正方法＊

方圣恩，張秋虎，林友勤，張笑華

（1.福州大學(xué)土木工程學(xué)院，福建 福州350108；2.合肥工業(yè)大學(xué)土木與水利工程學(xué)院，安徽 合肥230009）

引 言

過(guò)去20年來(lái)，有限元模型修正技術(shù)已經(jīng)廣泛應(yīng)用于參數(shù)識(shí)別和損傷識(shí)別等領(lǐng)域［1，2］。當(dāng)前絕大多數(shù)的模型修正方法都屬于確定性方法，無(wú)法考慮結(jié)構(gòu)參數(shù)和響應(yīng)的不確定性，使得實(shí)用性大大降低。由于模型誤差、構(gòu)件幾何尺寸和材料變異性、制造誤差、噪聲等因素的影響，不確定性廣泛存在于實(shí)際結(jié)構(gòu)的參數(shù)和響應(yīng)中［3］，是無(wú)法回避的問(wèn)題。因此，為了提高參數(shù)識(shí)別、損傷識(shí)別結(jié)果的可靠性，在有限元模型修正過(guò)程中結(jié)合概率方法進(jìn)行分析是十分必要的［4］。

工程分析中可以將不確定性分為偶然型不確定性（aleatory uncertainty）和認(rèn)知型不確定性（epistemic uncertainty）［5］。前者一般指結(jié)構(gòu)或構(gòu)件由于制造誤差或材料離散性所導(dǎo)致的幾何參數(shù)和材料變異性，可以被量化，但無(wú)法避免或人為消除；后者源于對(duì)結(jié)構(gòu)的認(rèn)識(shí)或測(cè)量信息上的不足（比如僅有少量的測(cè)量數(shù)據(jù)），但可以通過(guò)補(bǔ)充相關(guān)信息來(lái)降低甚至消除。在模型修正方面，通過(guò)結(jié)合概率分析方法可以有效考慮上述不確定性，使得修正結(jié)果的可靠性得到大大提高［6］。早期的概率模型修正方法主要針對(duì)測(cè)量誤差所引起的不確定性，通過(guò)最小化測(cè)量噪聲的方差來(lái)尋求參數(shù)均值的最佳估計(jì)值［6］。此外，貝葉斯方法的應(yīng)用可以獲得更加可靠的識(shí)別結(jié)果［7，8］，雖然修正過(guò)程會(huì)變得更為復(fù)雜、計(jì)算量更大。近些年來(lái)，采用攝動(dòng)方法在參數(shù)設(shè)計(jì)點(diǎn)上基于截?cái)嗵├照归_(kāi)的方式來(lái)表示修正方程的各個(gè)項(xiàng)，并以此構(gòu)建概率模型修正（probabilistic model updating）過(guò)程，可以有效提高修正效率［9］，因此得到了一定的應(yīng)用。但概率修正方法的精度依賴于對(duì)結(jié)構(gòu)參數(shù)和響應(yīng)概率分布特性的準(zhǔn)確估計(jì)，要求有大量的測(cè)試信息來(lái)建立準(zhǔn)確的概率分布函數(shù)。對(duì)大多數(shù)工程問(wèn)題而言，考慮到測(cè)試成本和可行性，上述要求往往無(wú)法得到滿足。此外，某些情況下參數(shù)或響應(yīng)的上下界（極限值）對(duì)工程分析或設(shè)計(jì)來(lái)說(shuō)更有實(shí)際意義，而非關(guān)注于具體的概率分布形式。此時(shí)可以采用區(qū)間數(shù)來(lái)表示結(jié)構(gòu)的參數(shù)和響應(yīng)，即構(gòu)建區(qū)間模型修正（interval model updating）問(wèn)題［10］。

目前為止，針對(duì)區(qū)間模型修正問(wèn)題的相關(guān)研究還非常少。同時(shí)由于區(qū)間數(shù)運(yùn)算法則和傳統(tǒng)數(shù)學(xué)算法非常不同［11］，因此主要還是在確定性框架內(nèi)進(jìn)行區(qū)間反向優(yōu)化（interval inverse optimization）問(wèn)題的求解，并不是真正意義上的區(qū)間優(yōu)化反演過(guò)程。對(duì)簡(jiǎn)單結(jié)構(gòu)而言，可以通過(guò)求解區(qū)間線性方程組來(lái)估計(jì)區(qū)間參數(shù)［12］。更常見(jiàn)的辦法是將區(qū)間模型修正問(wèn)題分解為兩個(gè)確定性修正過(guò)程，即分別獨(dú)立修正參數(shù)的上界和下界（或者區(qū)間中值和區(qū)間半徑），以避免區(qū)間數(shù)運(yùn)算可能導(dǎo)致的優(yōu)化困難［13］。此時(shí)可以采用一階泰勒展開(kāi)式來(lái)表示結(jié)構(gòu)的特征矩陣，并將區(qū)間模態(tài)參數(shù)作為目標(biāo)響應(yīng)。當(dāng)修正得到的響應(yīng)區(qū)間和實(shí)測(cè)區(qū)間基本吻合時(shí)，即可認(rèn)為優(yōu)化過(guò)程收斂，其間可以采用傳統(tǒng)的優(yōu)化算法進(jìn)行求解［14］。值得注意的是，區(qū)間算法的運(yùn)算過(guò)程普遍存在著一個(gè)問(wèn)題，即區(qū)間擴(kuò)張（intervals overestimation）現(xiàn)象。該問(wèn)題是由于區(qū)間變量在區(qū)間方程中多次出現(xiàn)所導(dǎo)致的，可以采用頂點(diǎn)法（vertex solution），通過(guò)合理組合參數(shù)的上下界求解響應(yīng)區(qū)間［15］。頂點(diǎn)法對(duì)求解簡(jiǎn)單結(jié)構(gòu)的特征值問(wèn)題（頻率）非常有效［16］，但其要求結(jié)構(gòu)的剛度和質(zhì)量矩陣是參數(shù)的線性函數(shù)。該條件對(duì)于復(fù)雜結(jié)構(gòu)而言通常無(wú)法滿足，因此是頂點(diǎn)法的重要應(yīng)用缺陷。與此同時(shí)，該方法也不能采用特征向量（模態(tài)振型）作為響應(yīng)，因?yàn)榇藭r(shí)無(wú)法保證參數(shù)－響應(yīng)關(guān)系的單調(diào)性［17］。因此，有必要尋求一種更為通用的方法求解區(qū)間模型修正問(wèn)題。此外，區(qū)間模型修正過(guò)程的計(jì)算效率也值得關(guān)注。為了提高修正效率，可以采用元模型（meta－model，如響應(yīng)面模型和Kriging模型）替代有限元模型，以建立結(jié)構(gòu)參數(shù)和響應(yīng)之間的關(guān)系式。然后通過(guò)不斷調(diào)整參數(shù)區(qū)間的大小，并在每個(gè)迭代步開(kāi)始時(shí)重構(gòu)Kriging模型實(shí)現(xiàn)修正過(guò)程［17］。這種方法具有比頂點(diǎn)法更優(yōu)的通用性，但其修正過(guò)程比較復(fù)雜且計(jì)算量大，不利于在復(fù)雜結(jié)構(gòu)問(wèn)題上的應(yīng)用。

本文將傳統(tǒng)響應(yīng)面模型和區(qū)間算法相結(jié)合，提出了區(qū)間響應(yīng)面模型（interval response surface model，縮寫IRSM）的概念，使得區(qū)間模型修正過(guò)程中可以直接在IRSM表達(dá)式上進(jìn)行區(qū)間數(shù)運(yùn)算和優(yōu)化計(jì)算，避免了將問(wèn)題分解為數(shù)個(gè)確定性過(guò)程來(lái)分別求解，從而大幅提高了修正效率。此外，本文方法對(duì)結(jié)構(gòu)參數(shù)－響應(yīng)之間的關(guān)系（單調(diào)或非單調(diào)，線性或非線性）沒(méi)有要求，且可以同時(shí)采用特征值和特征向量作為目標(biāo)響應(yīng)。更重要的是，采用區(qū)間響應(yīng)面模型不會(huì)發(fā)生區(qū)間擴(kuò)張現(xiàn)象，可以有效提高修正結(jié)果的可信度。最后，文中通過(guò)一個(gè)數(shù)值質(zhì)量－彈簧系統(tǒng)和一組試驗(yàn)鋼板驗(yàn)證了所提出方法的可行性和可靠性。

1 區(qū)間響應(yīng)面模型

1.1 區(qū)間擴(kuò)張問(wèn)題

區(qū)間數(shù)的運(yùn)算法則（包括加、減、乘、除）和傳統(tǒng)數(shù)學(xué)算法非常不同［11］，導(dǎo)致在區(qū)間反問(wèn)題的求解過(guò)程中直接進(jìn)行區(qū)間運(yùn)算非常困難，比如計(jì)算區(qū)間函數(shù)的梯度。因此，當(dāng)前絕大多數(shù)研究仍然在確定性框架下進(jìn)行求解，并非完全意義上的區(qū)間模型修正。此外，區(qū)間運(yùn)算過(guò)程中的區(qū)間擴(kuò)張現(xiàn)象始終存在，很大程度上降低了修正結(jié)果的可信度。舉例來(lái)說(shuō)，定義函數(shù)f（x），g（x）和h（x）如下

可以看到，3個(gè)函數(shù)從傳統(tǒng)數(shù)學(xué)理論上來(lái)說(shuō)是完全一樣的，僅表達(dá)式不同而已。若將函數(shù)中的變量x用區(qū)間數(shù)xI＝ ［0，1］表示，此時(shí)通過(guò)區(qū)間數(shù)運(yùn)算可得

可見(jiàn)3個(gè)函數(shù)求解得到的區(qū)間值完全不同，其中僅h（x）給出了正確解，而其他2個(gè)函數(shù)都出現(xiàn)了區(qū)間擴(kuò)張現(xiàn)象。這是由于x在f（x）和g（x）表達(dá)式中出現(xiàn)了2次，在運(yùn)算過(guò)程中被認(rèn)為是2個(gè)獨(dú)立變量所造成的。因此，文獻(xiàn)［15］提出了頂點(diǎn)法來(lái)解決這個(gè)問(wèn)題。但如前所述，頂點(diǎn)法在復(fù)雜問(wèn)題上的應(yīng)用有著很大局限性，因此有必要尋求更為通用的求解方法。

1.2 區(qū)間響應(yīng)面建模

本文將傳統(tǒng)響應(yīng)面模型和區(qū)間算法相結(jié)合，提出了IRSM的概念。傳統(tǒng)的響應(yīng)面模型實(shí)際上是一種顯式的多項(xiàng)式數(shù)學(xué)表達(dá)式，可以將系統(tǒng)的參數(shù)x和響應(yīng)y聯(lián)系起來(lái)［18］

式中β為常系數(shù)，通過(guò)對(duì)函數(shù)近似誤差ε的最小二乘估計(jì)來(lái)求解；k為參數(shù)個(gè)數(shù)；xi表示參數(shù)主效應(yīng)，xixj表示參數(shù)間相互效應(yīng)。式（3）為2階響應(yīng)面模型的完整表達(dá)式，其對(duì)大多數(shù)工程問(wèn)題來(lái)說(shuō)精度可以滿足要求，因而得到了廣泛應(yīng)用［19］。而在構(gòu)建IRSM過(guò)程中，需要對(duì)式（3）進(jìn)行改造。首先要忽略掉相互效應(yīng)項(xiàng)，因?yàn)榇蠖鄶?shù)情況下，參數(shù)間相互效應(yīng)對(duì)響應(yīng)面模型總方差的貢獻(xiàn)非常?。?0，21］，忽略后并不影響模型的精度，反而還能提高求解效率。然后利用配方法將式（3）轉(zhuǎn)化為完全平方項(xiàng)

這樣可以使表達(dá)式中xi僅出現(xiàn)一次，從而在區(qū)間運(yùn)算過(guò)程中避免了區(qū)間擴(kuò)張現(xiàn)象的出現(xiàn)。最后，將xi用區(qū)間數(shù)表示，可以得到

式中 上標(biāo)I表示區(qū)間數(shù)或區(qū)間變量。

式（5）即IRSM的理論表達(dá)式。利用該模型可以直接進(jìn)行區(qū)間運(yùn)算，以方便地求得響應(yīng)yI的區(qū)間值，并在優(yōu)化過(guò)程進(jìn)行快速反向求解。

1.3 區(qū)間響應(yīng)面和傳統(tǒng)響應(yīng)面的異同

本文提出的IRSM是建立在傳統(tǒng)的2階響應(yīng)面模型基礎(chǔ)上，通過(guò)對(duì)后者表達(dá)式的變換并引入?yún)^(qū)間數(shù)，得到IRSM表達(dá)式。該表達(dá)式的特點(diǎn)是每個(gè)區(qū)間參數(shù)都僅出現(xiàn)一次，以此有效避免區(qū)間運(yùn)算過(guò)程中的區(qū)間擴(kuò)張現(xiàn)象。此外，傳統(tǒng)響應(yīng)面僅能針對(duì)實(shí)數(shù)進(jìn)行分析，而IRSM可以同時(shí)考慮區(qū)間數(shù)和實(shí)數(shù)，即將實(shí)數(shù)作為“點(diǎn)區(qū)間”來(lái)進(jìn)行運(yùn)算，這也是IRSM的優(yōu)勢(shì)之一?？傮w來(lái)說(shuō)，這兩種響應(yīng)面類型之間既有聯(lián)系，也有區(qū)別，同時(shí)二者的應(yīng)用領(lǐng)域也有所不同。

2 區(qū)間模型修正過(guò)程

區(qū)間模型修正過(guò)程與傳統(tǒng)的確定性修正有著很大不同，前者需要采用區(qū)間算法來(lái)處理區(qū)間變量。除了前述的IRSM，本文提出的區(qū)間模型修正方法還需要考慮的問(wèn)題包括：1）參數(shù)和響應(yīng)的選??；2）區(qū)間目標(biāo)函數(shù)的建立；3）區(qū)間優(yōu)化算法。

2.1 參數(shù)和響應(yīng)的選取

選擇合適的參數(shù)和響應(yīng)是區(qū)間模型修正成功與否的關(guān)鍵一環(huán)。在結(jié)構(gòu)參數(shù)方面，一般可以選取代表結(jié)構(gòu)剛度和質(zhì)量的參數(shù)，比如彈性模量、截面慣性矩、材料密度等。在響應(yīng)的選擇上，可以采用結(jié)構(gòu)的靜動(dòng)力測(cè)量值如撓度、應(yīng)變、頻率、模態(tài)振型等作為目標(biāo)響應(yīng)，其中使用較多的是頻率（特征值）和振型（特征向量）。測(cè)量頻率要比振型簡(jiǎn)單得多，精度往往也更高。但振型包含了結(jié)構(gòu)的空間信息，對(duì)損傷定位來(lái)說(shuō)又是非常必要的，雖然考慮振型會(huì)使得區(qū)間模型修正問(wèn)題變得更為復(fù)雜。同時(shí)要注意的是，需要根據(jù)具體的修正策略和方法來(lái)確定響應(yīng)，比如頂點(diǎn)法就無(wú)法采用模態(tài)振型作為響應(yīng)，這點(diǎn)可以通過(guò)理論推導(dǎo)來(lái)說(shuō)明。

特征值對(duì)結(jié)構(gòu)質(zhì)量和剛度矩陣（元素）的偏微分可以表示為［22］：

式中 λ和φ分別表示特征值和特征向量；K和M表示結(jié)構(gòu)的剛度和質(zhì)量矩陣，其分塊矩陣（或元素）為kj和mj。因?yàn)镸和K分別為正定和半正定矩陣，所以式（6），（7）的符號(hào)在特定參數(shù)區(qū)間內(nèi)不會(huì)發(fā)生變化，即滿足單調(diào)性條件。因此，特征值區(qū)間的上下界容易通過(guò)剛度和質(zhì)量元素的上下界值組合求得。對(duì)特征值而言，其偏微分可以表示為［17］

其中

由式（9）可見(jiàn)，由于λi－λk的正負(fù)號(hào)會(huì)發(fā)生變化，即αik的正負(fù)號(hào)會(huì)發(fā)生變化，使得式（8）的單調(diào)性條件無(wú)法得到保證。因此，剛度和質(zhì)量元素的上下界值組合未必對(duì)應(yīng)特征向量區(qū)間的上下界，導(dǎo)致頂點(diǎn)法不適用。但是若采用區(qū)間響應(yīng)面模型進(jìn)行修正，則不要求滿足單調(diào)性條件，因此可以同時(shí)采用特征值和特征向量作為目標(biāo)響應(yīng)，這也是本文方法的主要優(yōu)點(diǎn)之一。

2.2 區(qū)間目標(biāo)函數(shù)

區(qū)間模型修正目標(biāo)函數(shù)的構(gòu)造與傳統(tǒng)的確定性修正問(wèn)題有所不同，前者不能直接采用區(qū)間響應(yīng)相對(duì)誤差的最小化來(lái)求解優(yōu)化問(wèn)題。舉例來(lái)說(shuō)，數(shù)值估計(jì)的區(qū)間響應(yīng)yI＝ ［20 ，30］ 與實(shí)測(cè)響應(yīng)區(qū)間＝［20，30］的相對(duì)誤差應(yīng)該為零。但對(duì)區(qū)間運(yùn)算來(lái)說(shuō)并非如此，比如：

由上式可見(jiàn)，在確定性模型修正中常用的目標(biāo)函數(shù)表達(dá)式F1和F2，經(jīng)過(guò)區(qū)間運(yùn)算后，其相對(duì)誤差值都不為零，因而不能直接用于區(qū)間模型修正中。有鑒于此，本文提出采用響應(yīng)yI的上下界和來(lái)構(gòu)建目標(biāo)函數(shù)，即有

式中m表示響應(yīng)的數(shù)目。

2.3 區(qū)間優(yōu)化算法

構(gòu)造目標(biāo)函數(shù)后，即可以通過(guò)最小化F求解區(qū)間優(yōu)化問(wèn)題

每個(gè)迭代步中首先調(diào)整參數(shù)區(qū)間的大??；然后利用區(qū)間響應(yīng)面模型快速計(jì)算結(jié)構(gòu)響應(yīng)區(qū)間，并與實(shí)測(cè)區(qū)間進(jìn)行對(duì)比，若不符合收斂性要求，則重新調(diào)整參數(shù)區(qū)間大小，再次優(yōu)化；反復(fù)進(jìn)行上述調(diào)整直至F值小于預(yù)先設(shè)置的容差（tolerance），此時(shí)迭代過(guò)程收斂。要注意的是，優(yōu)化過(guò)程要設(shè)置約束條件，即參數(shù)的下界x不能大于其上界，以保證修正參數(shù)區(qū)間的合理性。

在算法實(shí)現(xiàn)上，本文采用單目標(biāo)優(yōu)化方式，函數(shù)F中每階響應(yīng)在優(yōu)化過(guò)程中的權(quán)重相同，并通過(guò)建立線性約束來(lái)設(shè)置參數(shù)變化的界限，最后通過(guò)二次規(guī)劃算法求解F最小值的最優(yōu)解。值得一提的是，待修正參數(shù)的初始區(qū)間范圍（不確定性大?。┛梢栽O(shè)置較大，其對(duì)參數(shù)區(qū)間的最終估計(jì)精度影響很小。同時(shí)本文建議一次修正的參數(shù)數(shù)目不宜過(guò)多，因?yàn)轫憫?yīng)面構(gòu)建時(shí)所需的樣本數(shù)會(huì)隨著參數(shù)數(shù)目的增加而呈指數(shù)級(jí)增長(zhǎng)，使得計(jì)算量大為增加。同時(shí)過(guò)多參數(shù)也容易導(dǎo)致優(yōu)化過(guò)程出現(xiàn)病態(tài)問(wèn)題，使得優(yōu)化求解無(wú)法收斂或者降低了參數(shù)區(qū)間估計(jì)的精度。

圖1 3自由度質(zhì)量－彈簧系統(tǒng)Fig.1 Three degree of freedom mass－spring system

3 數(shù)值算例

本文首先通過(guò)一個(gè)數(shù)值質(zhì)量－彈簧系統(tǒng)來(lái)驗(yàn)證所提出的方法，如圖1所示。系統(tǒng)由3個(gè)質(zhì)量塊和6組彈簧組成，其中假定k1，k2，k5為不確定的區(qū)間參數(shù)。各參數(shù)的真實(shí)值或區(qū)間如下所設(shè)：

本例中采用系統(tǒng)的3階特征值λ1，λ2，λ3以及第一階特征向量的第一分量絕對(duì)值φ（1，1）作為目標(biāo)響應(yīng)，對(duì)應(yīng)地建立了4個(gè)IRSM，其中f1和φ（1，1）的表達(dá)式舉例如下：

可以看到，表達(dá)式中的變量均為區(qū)間參數(shù)，因此區(qū)間運(yùn)算可以直接在表達(dá)式上進(jìn)行，以快速計(jì)算響應(yīng)區(qū)間。同時(shí)可以基于上述模型，通過(guò)反向優(yōu)化過(guò)程來(lái)估計(jì)參數(shù)區(qū)間。本例中3個(gè)待修正參數(shù)的初始區(qū)間設(shè)為k1＝k2＝k5＝ ［0.5，1.5］，修正后的結(jié)果列于表1。由表1可見(jiàn)，區(qū)間修正后參數(shù)的平均誤差由初始的［37.5，25.0］%降低至［0.4，0.3］%，說(shuō)明參數(shù)區(qū)間的估計(jì)是非常精確的。圖2為系統(tǒng)特征值和特征向量的收斂圖，可以看到預(yù)測(cè)的響應(yīng)區(qū)間和實(shí)際區(qū)間吻合良好。由本算例分析結(jié)果可知，本文方法可以很好地估計(jì)不確定性參數(shù)的區(qū)間范圍，并且可以同時(shí)采用系統(tǒng)特征值和特征向量作為目標(biāo)響應(yīng)，比頂點(diǎn)法的應(yīng)用范圍更廣。此外，采用區(qū)間響應(yīng)面模型替代有限元模型，可以快速方便地進(jìn)行區(qū)間運(yùn)算，有效地簡(jiǎn)化了區(qū)間優(yōu)化問(wèn)題的復(fù)雜程度，并提高了修正效率。

表1 質(zhì)量－彈簧系統(tǒng)參數(shù)區(qū)間估計(jì)值Tab.1 Estimated parameter intervals of the mass－springsystem

圖2 質(zhì)量－彈簧系統(tǒng)響應(yīng)區(qū)間收斂圖Fig.2 Convergence plots of response intervals of the mass－spring system

4 試驗(yàn)驗(yàn)證

為了進(jìn)一步驗(yàn)證所提出的方法，本文實(shí)測(cè)了一組共55塊名義上完全相同的鋼板，如圖3所示。鋼板名義幾何尺寸為600mm（長(zhǎng)）×120mm（寬）×3 mm（厚）；名義材料特性為楊氏模量210GPa，剪切模量83GPa，質(zhì)量密度7 860kg／m3。試驗(yàn)中在鋼板一端的兩個(gè)角部（振幅最大處）各布置一個(gè)加速度傳感器，并采用橡皮筋懸掛的方式模擬自由邊界條件。然后對(duì)每塊鋼板進(jìn)行錘擊試驗(yàn)，同時(shí)測(cè)量錘擊信號(hào)和響應(yīng)信號(hào)，以獲取測(cè)點(diǎn)處的頻響函數(shù)，并提取鋼板前5階自振頻率，用于區(qū)間參數(shù)識(shí)別。試驗(yàn)過(guò)程為了盡可能減少人為因素和環(huán)境因素的干擾，所有鋼板都由同一個(gè)試驗(yàn)人員進(jìn)行測(cè)試，采用完全相同的試驗(yàn)步驟和試驗(yàn)儀器，并且保證鋼板的自由懸掛點(diǎn)和測(cè)試敲擊點(diǎn)也完全相同。同時(shí)，試驗(yàn)在相對(duì)恒溫的室內(nèi)進(jìn)行，以避免溫差造成的頻率變化。最后，由圖4可見(jiàn)，在最終獲取的鋼板模態(tài)中，第1，2，4階為彎曲模態(tài)，第3，5階為扭轉(zhuǎn)模態(tài)。

為了進(jìn)行區(qū)間模型修正，首先采用殼單元建立鋼板的有限元模型，沿鋼板短邊劃分10個(gè)單元，長(zhǎng)邊劃分30個(gè)單元，共得到300個(gè)殼單元，如圖4所示。然后根據(jù)鋼板的幾何尺寸和材料特性名義值計(jì)算得到前5階頻率的區(qū)間中值為44.37，122.89，136.54，242.11及279.80Hz。而對(duì)5階模態(tài)振型來(lái)說(shuō)，彎曲模態(tài)主要由鋼板厚度t和楊氏模量E控制；扭轉(zhuǎn)模態(tài)主要由t和剪切模量G控制?？紤]了扭轉(zhuǎn)模態(tài)后，區(qū)間參數(shù)識(shí)別的難度增加，更能有效驗(yàn)證所提出方法的可靠性。

圖3 鋼板錘擊試驗(yàn)Fig.3 Impact testing of the steel plates

圖4 鋼板有限元模型和振型Fig.4 The FE model and modal shapes of the plates

4.1 參數(shù)和響應(yīng)的選取

本算例選取t，E，G作為待修正參數(shù)，主要考慮到：1）鋼板在冷軋制造過(guò)程中，厚度t的變異性比長(zhǎng)度和寬度大得多，因?yàn)榍懈钿摪宓臅r(shí)候長(zhǎng)寬容易控制，而軋制過(guò)程則很難控制鋼板厚度（僅3mm）的均勻性（見(jiàn)表2所示）。更重要的是，同樣百分比的厚度改變所造成的鋼板頻率變化要比長(zhǎng)寬改變所造成的大得多。2）對(duì)于鋼材來(lái)說(shuō)，由于材料的均質(zhì)性，其材料特性是穩(wěn)定的。但加工工藝等因素可能會(huì)導(dǎo)致鋼板的材性呈現(xiàn)不均勻變化。

表2 厚度參數(shù)區(qū)間估計(jì)值Tab.2 Interval estimations of the thickness parameters

此外，考慮到頻率的測(cè)量精度相對(duì)較高，本例中選取鋼板的前5階頻率作為目標(biāo)響應(yīng)。由表3可見(jiàn)，對(duì)此類均質(zhì)鋼板來(lái)說(shuō)，頻率的不確定性很小，體現(xiàn)在實(shí)測(cè)頻率的小區(qū)間范圍上，因而對(duì)區(qū)間模型修正具有一定的難度。

表3 鋼板頻率區(qū)間修正值（厚度）Tab.3 Frequency intervals of the steel plates（thickness）

4.2 厚度t的區(qū)間識(shí)別

首先識(shí)別厚度t的區(qū)間范圍。為了提高識(shí)別精度，將t沿板長(zhǎng)邊方向分為3個(gè)等分區(qū)域t1，t2和t3。試驗(yàn)前對(duì)每個(gè)區(qū)域的厚度用高精度千分尺測(cè)量2次并取平均值，進(jìn)而得到55塊鋼板的厚度區(qū)間為t1＝ ［2.936，2.987］mm，t2＝ ［2.942，2.987］mm，t3＝ ［2.932，2.994］mm。

為了建立區(qū)間響應(yīng)面模型并進(jìn)行修正，t1，t2，t3的初始區(qū)間均設(shè)為［2.8，3.1］mm，然后對(duì)應(yīng)于5階頻率分別建立了5個(gè)IRSM。以第1階為例，其表達(dá)式如下

最后通過(guò)區(qū)間修正估計(jì)的參數(shù)區(qū)間值列于表2。由修正結(jié)果可見(jiàn)，厚度的估計(jì)值非常接近于實(shí)測(cè)值，絕對(duì)平均誤差由初始的［4.6，3.7］%降至［0.9，0.4］%，說(shuō)明區(qū)間優(yōu)化過(guò)程沒(méi)有發(fā)生區(qū)間擴(kuò)張現(xiàn)象，修正精度高。此外，修正后的頻率區(qū)間誤差（見(jiàn)表3所示）也由初始的［4.2，4.0］%降低至［0.7，0.8］%，說(shuō)明區(qū)間優(yōu)化過(guò)程是收斂的，如圖5所示。

圖5 鋼板頻率區(qū)間收斂圖（厚度參數(shù)）Fig.5 Convergence plots of frequency intervals of the steel plates（thickness parameters）

4.3 材料參數(shù)的區(qū)間識(shí)別

本算例同時(shí)對(duì)材料參數(shù)E，G的區(qū)間值進(jìn)行了識(shí)別，即認(rèn)為鋼板頻率的不確定性主要由材料參數(shù)的變異性所引起的。試驗(yàn)中不可能對(duì)所有鋼板進(jìn)行取樣，因此僅對(duì)幾塊鋼板材料試樣進(jìn)行了單向拉伸試驗(yàn)，并根據(jù)試驗(yàn)結(jié)果設(shè)置了E，G的初始區(qū)間分別為E＝ ［190，220］GPa，G＝［77，89］GPa。隨后建立了前5階頻率的IRSM，其中第1階的表達(dá)式舉例如下

經(jīng)過(guò)區(qū)間修正得到E，G的區(qū)間估計(jì)值分別為E＝［196.5，203.6］GPa，G＝［79.5，83.4］GPa，比初始區(qū)間大為縮小。二者的區(qū)間中值分別為200.1 GPa和81.5GPa，屬于鋼材材料特性的合理值，說(shuō)明參數(shù)識(shí)別結(jié)果是合理的。修正后的頻率區(qū)間誤差（如表4所示）也由初始的［1.9，4.1］%降低至［0.4，0.2］%，說(shuō)明區(qū)間優(yōu)化過(guò)程是收斂的，如圖6所示。

表4 鋼板頻率區(qū)間修正值（材料參數(shù)）Tab.4 Frequency intervals of the steel plates（material properties）

圖6 鋼板頻率區(qū)間收斂圖（材料參數(shù)）Fig.6 Convergence plots of frequency intervals of the steel plates（material parameters）

5 結(jié) 論

目前在模型修正領(lǐng)域，針對(duì)不確定性參數(shù)的非概率識(shí)別方法還很少，特別是基于區(qū)間分析的相關(guān)研究。本文將響應(yīng)面模型和區(qū)間算法相結(jié)合，提出了區(qū)間響應(yīng)面模型的概念，用于簡(jiǎn)化區(qū)間模型修正過(guò)程，以避免區(qū)間擴(kuò)張問(wèn)題并提高修正效率。所提出的方法可以精確識(shí)別一個(gè)數(shù)值質(zhì)量－彈簧系統(tǒng)的彈簧剛度區(qū)間范圍，同時(shí)也能很好地識(shí)別出一組試驗(yàn)鋼板的幾何及材料參數(shù)的變異性區(qū)間。本文的研究可以得到以下結(jié)論：

1）利用IRSM可以快速、方便地估計(jì)結(jié)構(gòu)響應(yīng)的區(qū)間，同時(shí)還有利于區(qū)間優(yōu)化反問(wèn)題的求解，避免了區(qū)間擴(kuò)張現(xiàn)象和參數(shù)靈敏度的計(jì)算，大幅簡(jiǎn)化了優(yōu)化過(guò)程并提高了修正效率。

2）區(qū)間響應(yīng)面模型修正過(guò)程可以同時(shí)采用特征值和特征向量作為目標(biāo)響應(yīng)，不要求參數(shù)－響應(yīng)的單調(diào)性關(guān)系，因此比頂點(diǎn)法的應(yīng)用范圍更廣泛。

3）區(qū)間響應(yīng)面模型修正所需的目標(biāo)函數(shù)可以根據(jù)響應(yīng)的上下界建立，而非直接采用區(qū)間數(shù)，以此避免區(qū)間運(yùn)算導(dǎo)致的誤差。

4）所提出的方法為區(qū)間模型修正技術(shù)提供了一個(gè)新的問(wèn)題解決途徑。

［1］ Friswell M I，Mottershead J E.Finite Element Model Updating in Structural Dynamics［M］.Kluwer Academic Publishers，Dordrecht，Netherlands，1995.

［2］ Fang S E，Perera R，Roeck G D.Damage identification of a reinforced concrete frame by finite element model updating using damage parameterization ［J］.Journal of Sound and Vibration，2008，313：544—559.

［3］ Alvin K，Oberkampf W，Diegert K，et al.Uncertainty quantification in computational structural dynamics：a new paradigm for model validation［A］.Proceedings of the 16th International Modal Analysis Conference［C］，1998：1 191—1 197.

［4］ 馮新，李國(guó)強(qiáng)，周晶.土木工程結(jié)構(gòu)健康診斷中的統(tǒng)計(jì)識(shí)別方法綜述［J］.地震工程與工程振動(dòng)，2005，25（2）：105—113.Feng Xin，Li Guoqiang，Zhou Jing.State－of－the－art of statistical identification for structural health diagnosis in civil engineering ［J］.Earthquake Engineering and Engineering Vibration，2005，25（2）：105—113.

［5］ Moens D.A non－probabilistic finite element approach for structural dynamic analysis with uncertain parameters［D］.Catholic University of Leuven，Belgium，2002.

［6］ Friswell M I.The adjustment of structural parameters using a minimum variance estimator［J］.Mechanical Systems and Signal Processing，1989，3 （2）：143—155.

［7］ Beck J L，Katafygiotis L S.Updating models and their uncertainties.I：Bayesian statistical framework ［J］.Journal of Engineering Mechanics.1998，124 （4）：455—461.

［8］ 韓芳，鐘冬望，汪君.基于貝葉斯法的復(fù)雜有限元模型修正研究［J］.振動(dòng)與沖擊，2012，31（1）：39—43.HAN Fang，ZHONG Dongwang，WANG Jun.Complicated finite element model updating based on Bayesian method ［J］.Journal of Vibration and Shock，2012，31（1）：39—43.

［9］ Khodaparast H H，Mottershead J E，F(xiàn)riswell M I.Perturbation methods for the estimation of parameter variability in stochastic model updating［J］.Mechanical Systems and Signal Processing，2008，22（8）：1 751—1 773.

［10］王登剛，秦仙蓉.結(jié)構(gòu)計(jì)算模型修正的區(qū)間反演方法［J］.振動(dòng)工程學(xué)報(bào)，2004，17（2）：205—209.Wang Denggang，Qin Xianrong.Interval method for computational model updating of dynamic structures［J］.Journal of Vibration Engineering，2004，17（2）：205—209.

［11］Moore R，Kearfott R，Cloud M.Introduction to Interval Analysis［M］.Society for Industrial and Applied Mathematics，Philadelphia，2009.

［12］Rao S S，Berke L.Analysis of uncertain structural systems using interval analysis ［J］.AIAA Journal，1997 35（4）：727—735.

［13］Li S L，Li H，Ou J P.Model updating for uncertain structures with interval parameters［A］.Asia－Pacific Workshop on Structural Health Monitoring［C］.Yokohama，Japan，2006.

［14］Jiang C，Han X，Liu G.R.Optimization of structures with uncertain constraints based on convex model and satisfaction degree of interval［J］.Computer Methods in Applied Mechanics and Engineering，2007，196：4 791—4 800.

［15］Dong W，Shah H.Vertex method for computing functions of fuzzy variables［J］.Fuzzy Sets and Systems，1987，24（1）：65—78.

［16］Qiu Z P，Wang X J，F(xiàn)riswell M I.Eigenvalue bounds of structures with uncertain－but－bounded parameters［J］.Journal of Sound and Vibration，2005，282（1）：297—312.

［17］Khodaparast H H，Mottershead J E，Badcock K J.Interval model updating with irreducible uncertainty using the Kriging predictor ［J］.Mechanical Systems and Signal Processing，2011，25（4）：1 204—1 226.

［18］Myers R H，Montgomery D C，Anderson－Cook C M.Response Surface Methodology：Process and Product Optimization Using Designed Experiments［M］.3rd Edition，New York：John Wiley ＆Sons，2009.

［19］Cundy A L.Use of response surface metamodels in damage identification of dynamic structures［D］.Virginia Polytechnic Institute and State University，2002.

［20］Fang S E，Perera R.A response surface methodology based damage identification technique［J］.Smart Materials and Structures，2009，18（6）：065009.

［21］Fang S E，Perera R.Damage identification by response surface model updating using D－optimal design［J］.Mechanical Systems and Signal Processing，2011，25（2）：717—733.

［22］Fox R，Kapoor M.Rates of change of eigenvalues and eigenvectors ［J］.AIAA Journal，1968，6（12）：2 426—2 429.