楊秋偉， 周 聰， 李翠紅， 羅 帥

(1.紹興文理學(xué)院 土木工程系，紹興 312000； 2.寧波工程學(xué)院 建筑與交通工程學(xué)院，寧波 315211；3.寧波工程學(xué)院 浙江省土木工程工業(yè)化建造工程技術(shù)研究中心,寧波 315211)

1 引 言

有限單元法已廣泛用于土木、機(jī)械、航空航天、汽車和船舶等眾多工程設(shè)計(jì)領(lǐng)域中，成為解決復(fù)雜工程分析計(jì)算問(wèn)題的有效途徑。然而，由于大型工程結(jié)構(gòu)的復(fù)雜性，用有限元軟件建立的結(jié)構(gòu)模型與真實(shí)結(jié)構(gòu)總是存在一定的差異，具體體現(xiàn)在，由有限元模型計(jì)算所得的結(jié)構(gòu)響應(yīng)參數(shù)(如位移、頻率和振型等)總是與儀器測(cè)試所得的結(jié)構(gòu)響應(yīng)參數(shù)存在差異。當(dāng)這種差異較大時(shí)，則必須對(duì)結(jié)構(gòu)的有限元模型進(jìn)行修正，確保計(jì)算所得與測(cè)試所得的響應(yīng)參數(shù)盡可能接近。只有修正后的有限元模型方可用于力學(xué)分析與計(jì)算。進(jìn)一步，如果結(jié)構(gòu)中存在損傷，通過(guò)修正完好結(jié)構(gòu)的有限元模型來(lái)與測(cè)試所得的損傷結(jié)構(gòu)的響應(yīng)參數(shù)相匹配，根據(jù)模型的修正量就可以判斷出結(jié)構(gòu)的損傷單元及損傷程度，這就是基于模型修正的結(jié)構(gòu)損傷識(shí)別方法基本原理。近年來(lái)，模型修正法和以它為基礎(chǔ)的損傷識(shí)別方法已成為眾多工程技術(shù)領(lǐng)域的一項(xiàng)重要課題[1-4]。

殘余力向量法[5-16]是一類常用的損傷識(shí)別算法，由殘余力向量可以判定損傷位置，求解殘余力方程可以得到損傷程度。求解方法有最小范數(shù)法[5,6]、最小秩法[7,8]或者其他的迭代優(yōu)化算法[14-18]。目前，已有的殘余力向量法都是基于動(dòng)力測(cè)試的模態(tài)參數(shù)，因此可以統(tǒng)稱為動(dòng)力殘余力向量法。和動(dòng)力測(cè)試數(shù)據(jù)相比，靜力測(cè)試數(shù)據(jù)往往精度更高，且無(wú)需模態(tài)分析等復(fù)雜操作，易于實(shí)現(xiàn)。因此，本文基于靜力測(cè)試數(shù)據(jù)(位移參數(shù))，結(jié)合有限元?jiǎng)偠染仃?，定義了靜力殘余力向量，根據(jù)其中不為零的元素來(lái)判斷損傷自由度，再根據(jù)自由度和單元之間的對(duì)應(yīng)關(guān)系來(lái)確定損傷位置，在此基礎(chǔ)上，提出一種求解損傷參數(shù)的代數(shù)解法。另外，鑒于梁、剛架或?qū)嶓w結(jié)構(gòu)中節(jié)點(diǎn)角位移難以測(cè)量的情況，進(jìn)一步提出了一種靜力縮聚殘余力向量法，只需要應(yīng)用節(jié)點(diǎn)線位移數(shù)據(jù)，從而有效拓寬了所提方法的應(yīng)用范圍。數(shù)值算例結(jié)果說(shuō)明所提方法合理有效。

2 靜力殘余力向量法

一般而言，結(jié)構(gòu)的損傷主要導(dǎo)致結(jié)構(gòu)剛度的變化，質(zhì)量的變化基本可以忽略不計(jì)。不失一般性，假設(shè)K0和K分別為結(jié)構(gòu)損傷前后的剛度矩陣(n×n維)，由于損傷所導(dǎo)致的剛度損失量為ΔK，即

(1)

式中Ki為第i個(gè)單元?jiǎng)偠染仃?，αi為第i個(gè)單元的損傷參數(shù)(αi∈[0,1])，αi=0表示第i個(gè)單元未損傷，αi=1表示第i個(gè)單元完全損傷，N為單元體總數(shù)目。結(jié)構(gòu)損傷前后，在荷載向量l作用下的位移向量分別為

K0d0=l, (K0-ΔK)d=l

(2,3)

式中d0和d分別為完好結(jié)構(gòu)和損傷結(jié)構(gòu)的位移向量，其可以通過(guò)靜力測(cè)試得到。由方程(3)變形可得

ΔK·d={δ}， {δ}=K0d-l

(4，5)

方程(4)即為靜力殘余力方程，其中向量{δ}稱為靜力殘余力向量，由方程(5)可知，只要測(cè)試獲得損傷結(jié)構(gòu)在已知荷載向量l作用下的位移d，即可由方程(5)計(jì)算得出其靜力殘余力向量。根據(jù)方程(4)，由靜力殘余力向量{δ}來(lái)進(jìn)行損傷定位的原理在于，結(jié)構(gòu)損傷一般只導(dǎo)致ΔK的少數(shù)行向量不為零，即和損傷對(duì)應(yīng)的自由度所在的行向量不為零，因此，靜力殘余力向量{δ}中不為零的元素即對(duì)應(yīng)著損傷自由度，由損傷自由度并根據(jù)幾何關(guān)系即可確定發(fā)生損傷的具體單元。具體而言，可將方程(4)改寫為以下方程組，

(6)

當(dāng)損傷位置由上述方法判定出來(lái)以后，結(jié)合方程(1，6)，可以進(jìn)一步求解損傷程度，即求解損傷參數(shù)αi。分兩種情況討論。

首先討論單個(gè)損傷情況，如果由靜力殘余力向量{δ}判定出來(lái)只有一個(gè)單元體發(fā)生損傷，則只有一個(gè)未知的損傷參數(shù)αi需要求解，此時(shí)方程(1)簡(jiǎn)化為

ΔK=αiKi

(7)

相應(yīng)地，方程(6)簡(jiǎn)化為

(8)

(9)

方程(9)中δj取{δ}的任意一個(gè)非零元素即可。

(10)

方程(10)代入式(6)整理可得

Π{αi}′={δ}′

(11)

式中Π為整理所得的系數(shù)矩陣，{αi}′為q個(gè)未知損傷參數(shù)組成的列向量，{δ}′為q個(gè)非零的δj組成的列向量。由方程(11)即可求得損傷程度為

{αi}′=Π-1{δ}′

(12)

對(duì)于梁、剛架和實(shí)體結(jié)構(gòu)，其節(jié)點(diǎn)位移既有線位移也有角位移，實(shí)踐中由于角位移難以測(cè)量，上述的靜力殘余力向量法需要進(jìn)一步改良方可使用。為此，本文采用模型縮聚來(lái)對(duì)這些角位移進(jìn)行縮減處理，從而使得縮聚后的相關(guān)公式都只包含各節(jié)點(diǎn)線位移，拓寬了靜力殘余力向量法的應(yīng)用范圍?？s聚方法具體過(guò)程如下。

首先,將結(jié)構(gòu)的全部自由度分為平動(dòng)自由度(保留自由度)和轉(zhuǎn)動(dòng)自由度(縮減自由度)兩類，相應(yīng)地將方程(2)改寫為

(13)

式中上標(biāo)m和s分別表示平動(dòng)自由度和轉(zhuǎn)動(dòng)自由度，m和s同時(shí)也表示平動(dòng)自由度和轉(zhuǎn)動(dòng)自由度的數(shù)目。若所有的轉(zhuǎn)動(dòng)自由度上均未施加外荷載，即ls=0，則由方程(13)的第二行可得

(14)

由方程(14)可得

(15)

(16，17)

式中矩陣Θ為自由度轉(zhuǎn)換矩陣。將方程(16)代入式(13)，并兩邊左乘ΘT可得

Krdm=lr,Kr=ΘTK0Θ,lr=ΘTl

(18～20)

利用方程(18)縮聚后的剛度矩陣和荷載向量，按照前述靜力殘余力向量同樣的推導(dǎo)過(guò)程，可推導(dǎo)得到縮聚后的殘余力向量計(jì)算公式如下，

{δ}r=Krdm-lr

(21)

式中Kr為完好結(jié)構(gòu)縮聚后的剛度矩陣，lr為縮聚后的荷載向量，dm為損傷后結(jié)構(gòu)各平動(dòng)自由度處的線位移向量。相應(yīng)地，也可以得到模型縮聚后類似于方程(9，12)的計(jì)算損傷程度公式為

(22)

(23)

3 數(shù)值算例

3.1 桁架結(jié)構(gòu)

圖1 桁架結(jié)構(gòu)及靜力加載

圖2 單元10損傷時(shí)的靜力殘余力向量(單位：kN)

表1 損傷結(jié)構(gòu)位移和殘余力向量計(jì)算結(jié)果

接下來(lái)討論多個(gè)損傷情況，不失一般性，假設(shè)單元16和31剛度分別損傷15%和30%，這種情況下?lián)p傷結(jié)構(gòu)的位移數(shù)據(jù)和殘余力向量分別列于表1的第4列和第5列，為了更直觀，將表1的第5列數(shù)據(jù)繪入圖3。可以看出，殘余力不為零的自由度編號(hào)為9/13，25/26/31/32，根據(jù)圖1，自由度9和13對(duì)應(yīng)的節(jié)點(diǎn)為第5和第7號(hào)節(jié)點(diǎn)，顯然，這兩個(gè)節(jié)點(diǎn)之間的桿件剛好為單元16，而自由度25/26/31/32分別對(duì)應(yīng)著13和16號(hào)節(jié)點(diǎn)，這兩個(gè)節(jié)點(diǎn)之間為單元31。因此，根據(jù)殘余力向量中不為零的元素最終可以判定單元16和31為發(fā)生損傷的單元。接下來(lái)計(jì)算損傷程度，由于欲求的損傷參數(shù)個(gè)數(shù)為2，需要列兩個(gè)方程(方程(11))，這里選擇第9個(gè)自由度(對(duì)應(yīng)損傷單元16)和第32個(gè)自由度(對(duì)應(yīng)損傷單元31)來(lái)列方程，方程(11)的具體計(jì)算結(jié)果為

圖3 單元16和31損傷時(shí)的靜力殘余力向量(單位：kN)

(24)

進(jìn)一步由方程(24)計(jì)算可得損傷參數(shù)為α16=0.15和α31=0.3，剛好與假設(shè)值15%和30%相等。

需要說(shuō)明的是，在上述模擬過(guò)程中，并未考慮有限元模型誤差與測(cè)試數(shù)據(jù)誤差，因此計(jì)算所得的損傷程度值和假設(shè)值完全相等，這也直觀地說(shuō)明了本文所提的損傷程度代數(shù)解是一種理論上的精確解。如果考慮數(shù)據(jù)誤差，如對(duì)于單元10損傷情況，在損傷結(jié)構(gòu)各個(gè)位移數(shù)據(jù)中添加誤差水平為3%的隨機(jī)數(shù)，則可得靜力殘余力向量如圖4所示，可以看出第7和第11自由度(對(duì)應(yīng)單元10)處殘余力明顯偏大，故仍然可以判斷單元10為損傷單元，其損傷程度計(jì)算值為0.1989，與假設(shè)值0.2很接近。如果有些情況難以根據(jù)殘余力向量圖形做出準(zhǔn)確的損傷定位，可以結(jié)合殘余力向量圖形和損傷程度計(jì)算結(jié)果來(lái)綜合進(jìn)行判定，具體見接下來(lái)的梁結(jié)構(gòu)算例。

圖4 利用含噪聲數(shù)據(jù)計(jì)算所得靜力殘余力向量(單元10損傷，單位：kN)

3.2 梁結(jié)構(gòu)

圖5 梁結(jié)構(gòu)

不失一般性，假設(shè)單元15剛度損傷20%，圖6給出了無(wú)噪聲時(shí)靜力縮聚殘余力向量的計(jì)算結(jié)果?？梢钥闯?，絕對(duì)值偏大的殘余力集中于節(jié)點(diǎn)13，14，15和16處。根據(jù)單元體和平動(dòng)自由度編號(hào)之間的對(duì)應(yīng)關(guān)系，位于這4個(gè)節(jié)點(diǎn)之間的是第14，15和16號(hào)單元。由此可見，采用靜力縮聚殘余力向量，由于模型縮聚所帶來(lái)的誤差，不能精確定位于真正發(fā)生損傷的第15號(hào)單元處，而是定位于以損傷單元為中心的局部小區(qū)域。為了進(jìn)一步準(zhǔn)確判斷，分別計(jì)算第14，15和16號(hào)單元的損傷程度值。首先，計(jì)算第14號(hào)單元，該單元對(duì)應(yīng)的節(jié)點(diǎn)編號(hào)為13和14，這里采用節(jié)點(diǎn)13對(duì)應(yīng)的殘余力來(lái)進(jìn)行計(jì)算(由于節(jié)點(diǎn)14的殘余力可能由于第14號(hào)單元或第15號(hào)單元引起，故此處不宜采用節(jié)點(diǎn)14)，所得結(jié)果為α14=-0.3581，據(jù)此判斷第14號(hào)單元并非真正損傷單元。接下來(lái)，計(jì)算第16號(hào)單元，該單元對(duì)應(yīng)的節(jié)點(diǎn)編號(hào)為15和16，這里采用節(jié)點(diǎn)16對(duì)應(yīng)的殘余力進(jìn)行計(jì)算，結(jié)果為α16=-0.3402，據(jù)此判斷第16號(hào)單元也并非真正損傷單元。最后，計(jì)算第15號(hào)單元，該單元對(duì)應(yīng)的節(jié)點(diǎn)編號(hào)為14和15，為了可靠起見，分別采用這兩個(gè)節(jié)點(diǎn)處殘余力對(duì)應(yīng)的方程來(lái)進(jìn)行計(jì)算，所得結(jié)果分別為α15=0.2156(節(jié)點(diǎn)14)和α15=0.2159(節(jié)點(diǎn)15)，這兩個(gè)結(jié)果均說(shuō)明第15號(hào)單元為真正損傷單元，損傷程度計(jì)算值可取為二者的平均值，顯然與假設(shè)值0.2很接近?？紤]在位移數(shù)據(jù)中添加3%噪聲水平的隨機(jī)數(shù)，所得靜力縮聚殘余力向量如圖7所示，可見絕對(duì)值偏大的殘余力仍然集中于節(jié)點(diǎn)13，14，15和16處，利用和上述相同的損傷程度計(jì)算過(guò)程，所得結(jié)果分別為α14=-0.3005，α16=-0.3714，α15=0.2012(節(jié)點(diǎn)14)和α15=0.2235(節(jié)點(diǎn)15)，故仍然可以判定出第15號(hào)單元為損傷單元，說(shuō)明所提的靜力縮聚殘余力向量法是合理可行的。

圖6 無(wú)噪聲時(shí)靜力縮聚殘余力向量(單元15損傷，單位：kN)

圖7 有噪聲時(shí)靜力縮聚殘余力向量(單元15損傷，單位：kN)

4 結(jié) 論

提出一種靜力殘余力向量法用于結(jié)構(gòu)損傷評(píng)估，利用結(jié)構(gòu)靜力測(cè)試數(shù)據(jù)和有限元模型剛度矩陣，根據(jù)靜力殘余力向量中不為零的元素來(lái)判斷損傷自由度；再根據(jù)自由度和單元之間的對(duì)應(yīng)關(guān)系來(lái)確定發(fā)生損傷的位置，并提出一種求解損傷參數(shù)的代數(shù)解法；另外，針對(duì)實(shí)踐中角位移難以測(cè)量的情況，進(jìn)一步提出了一種靜力縮聚殘余力向量法，有效拓寬了所提方法的應(yīng)用范圍。數(shù)值算例結(jié)果表明，所提方法合理可行，對(duì)單個(gè)損傷和多個(gè)損傷情況均可適用，可為解決結(jié)構(gòu)損傷識(shí)別問(wèn)題提供參考。