秦世強(qiáng)， 甘耀威， 康俊濤

(武漢理工大學(xué) 土木工程與建筑學(xué)院， 武漢 430070)

建立一個(gè)準(zhǔn)確、能夠代表實(shí)際結(jié)構(gòu)行為的有限元模型對(duì)橋梁狀態(tài)評(píng)估和監(jiān)測(cè)有著重要的作用。然而，按照設(shè)計(jì)圖紙建立的初始有限元模型，由于模型簡(jiǎn)化、結(jié)構(gòu)設(shè)計(jì)參數(shù)的不確定性、結(jié)構(gòu)退化或損傷等因素的影響[1]，往往難以代表實(shí)際結(jié)構(gòu)行為。因此，需要對(duì)初始有限元模型進(jìn)行修正，以降低結(jié)構(gòu)響應(yīng)的有限元計(jì)算值與實(shí)測(cè)值之間的相對(duì)誤差。有限元模型修正的本質(zhì)是一個(gè)優(yōu)化問題，其目標(biāo)函數(shù)是由結(jié)構(gòu)響應(yīng)有限元計(jì)算值與實(shí)測(cè)值的殘差函數(shù)構(gòu)成。對(duì)于橋梁結(jié)構(gòu)模型修正的目標(biāo)函數(shù)，通常具有多維、高度非線性和多個(gè)局部極值點(diǎn)的特征；傳統(tǒng)搜索算法如梯度算法、牛頓迭代法在處理這類問題時(shí)，容易陷入局部最優(yōu)，難以獲得目標(biāo)函數(shù)的全局最優(yōu)解[2]。近年來基于智能優(yōu)化理論的粒子群算法[3]、遺傳算法[4]、猴群算法[5]、萬有引力搜索算法(gravitational search algorithm ，GSA)[6]等迅速發(fā)展，相比傳統(tǒng)搜索算法，智能進(jìn)化算法在處理非線性、多極值等問題時(shí)具有一定的優(yōu)勢(shì)。Alkayem等[7]對(duì)基于智能算法的有限元模型修正進(jìn)行了全面的綜述，總結(jié)了各類單目標(biāo)優(yōu)化算法和多目標(biāo)優(yōu)化算法在有限元模型修正中的應(yīng)用現(xiàn)狀；同時(shí)，作者指出提升智能算法的搜索精度和效率，仍是模型修正和損傷識(shí)別領(lǐng)域一個(gè)值得研究的方向。

GSA[8]是Esmat Rashedi等提出的一種元啟發(fā)式算法，是一種模擬物理學(xué)中萬有引力的新的優(yōu)化算法。它通過種群中粒子間的相互作用來指導(dǎo)群體進(jìn)行智能優(yōu)化搜索；已有研究表明[9]：GSA的收斂性明顯優(yōu)于粒子群算法、遺傳算法等仿生優(yōu)化算法。但是GSA依然存在著容易早熟、后期迭代速度慢、不易跳出局部最優(yōu)解等問題。針對(duì)GSA的上述問題，已有文獻(xiàn)提出了一些改進(jìn)措施。谷文祥等[10]提出一種新的局部搜索算法并結(jié)合GSA求解流水線調(diào)度問題，有效的避免算法陷入局部最優(yōu)值，提高了解的質(zhì)量，驗(yàn)證了所得算法的有效性；李超順等[11]結(jié)合粒子群算法的運(yùn)動(dòng)特點(diǎn)提出了改進(jìn)的引力搜索算法，并將其應(yīng)用到勵(lì)磁控制系統(tǒng)PID(proportional-integral-derivative)參數(shù)優(yōu)化，并與傳統(tǒng)群體優(yōu)化算法進(jìn)行了充分對(duì)比試驗(yàn)，驗(yàn)證了其優(yōu)化算法的有效性；李沛等[12]在GSA的速度更新部分引入粒子群算法中的記憶和群體信息交流功能，并提出了基于權(quán)值的粒子慣性質(zhì)量更新公式，驗(yàn)證了所提算法能有效規(guī)劃出無人機(jī)的最優(yōu)航路；蔣建國(guó)等[13]針對(duì)基本GSA在處理優(yōu)化問題時(shí)會(huì)出現(xiàn)發(fā)散的情況,通過限制粒子的速度同時(shí)更改算法中的參數(shù)，并將改進(jìn)的GSA運(yùn)用到邊坡穩(wěn)定性分析中，搜索出臨界滑動(dòng)面并結(jié)合極限平衡法計(jì)算出相應(yīng)的最小安全系數(shù)，并驗(yàn)證了其具有更好的優(yōu)化性能。本文在上述研究的基礎(chǔ)上，結(jié)合遺傳算法提出一種改進(jìn)的GSA，利用荷載試驗(yàn)結(jié)果并結(jié)合Kriging代理模型，實(shí)現(xiàn)了一座鋼混疊合梁組合式系桿拱體系橋梁的模型修正。

1 基于Kriging模型的模型修正

有限元模型修正的本質(zhì)是通過改變模型中的結(jié)構(gòu)設(shè)計(jì)參數(shù)，減小有限元計(jì)算值與實(shí)測(cè)值之間的誤差。因此，模型修正的目標(biāo)函數(shù)可以定義為

(1)

式中：x為待修正的設(shè)計(jì)參數(shù)；ri(x)為第i類結(jié)構(gòu)響應(yīng)計(jì)算值與實(shí)測(cè)值之間的殘差函數(shù)；m為模型修正中考慮的結(jié)構(gòu)響應(yīng)的類型種數(shù)；αi為不同殘差函數(shù)的權(quán)重系數(shù)。

為獲得式(1)所表達(dá)的目標(biāo)函數(shù)在設(shè)計(jì)空間內(nèi)的最小值，一般需通過迭代法求解。然而，由于結(jié)構(gòu)響應(yīng)與設(shè)計(jì)參數(shù)間隱函數(shù)的關(guān)系，使得每次迭代均需要調(diào)用有限元模型來計(jì)算殘差函數(shù)，這一過程十分耗時(shí)，特別是對(duì)復(fù)雜的橋梁結(jié)構(gòu)有限元模型而言?；诖砟Ｐ偷哪Ｐ托拚椒ㄒ虼说靡园l(fā)展；常用的代理模型包括多項(xiàng)式響應(yīng)面[14]、徑向基函數(shù)[15]、Kriging模型[16]等。有研究表明[17]：Kriging模型不僅可以給出設(shè)計(jì)空間內(nèi)任意一點(diǎn)的預(yù)測(cè)值，同時(shí)能給出該點(diǎn)的預(yù)測(cè)誤差?；贙riging模型的模型修正的主要流程可分為3個(gè)步驟:① 需要用初始有限元模型得到結(jié)構(gòu)的響應(yīng)和設(shè)計(jì)參數(shù)樣本，可通過試驗(yàn)設(shè)計(jì)[18]完成，常用的試驗(yàn)設(shè)計(jì)有中心復(fù)合設(shè)計(jì)、D最優(yōu)設(shè)計(jì)、均勻設(shè)計(jì)和拉丁超立方設(shè)計(jì)等；② 利用得到的響應(yīng)和設(shè)計(jì)參數(shù)樣本生成Kriging模型并做精度測(cè)試；③ 基于優(yōu)化算法尋找目標(biāo)函數(shù)在設(shè)計(jì)空間內(nèi)的最小值。

Kriging代理模型是由確定性的回歸模型和隨機(jī)過程兩部分組成，結(jié)構(gòu)響應(yīng)矩陣Y(x)可以用Kriging模型表達(dá)成設(shè)計(jì)參數(shù)x的表達(dá)式，公式為

(2)

2 基本萬有引力搜索算法

GSA是一種模擬萬有引力定律的智能搜索算法。GSA中每個(gè)粒子都可以看作一個(gè)個(gè)體，每個(gè)粒子都具有自身的位置和移動(dòng)速度，并且可以保留迄今為止粒子所能找到的最優(yōu)位置以及所有粒子當(dāng)前能找到的全局最優(yōu)位置。在標(biāo)準(zhǔn)的GSA中，粒子i和粒子j之間的引力可表示為

(3)

式中：pd,i(t)為在t時(shí)刻第i個(gè)個(gè)體在d維空間中所在的位置；Mpi(t)為t時(shí)刻受力個(gè)體i的質(zhì)量；Maj(t)為t時(shí)刻施力個(gè)體的質(zhì)量；ε為一個(gè)很小的常數(shù)；Rij(t)為粒子i與j之間的歐幾里得距離，可表示為

(4)

G(t)為在第t時(shí)刻的萬有引力常系數(shù)，可表示為

(5)

式中：G0為常系數(shù)的初始值；α為下降系數(shù)；T為總迭代次數(shù)。因此，第i個(gè)個(gè)體在d維空間中的合力可表示為

(6)

式中：N為粒子種群的個(gè)數(shù)；rand則為0到1的隨機(jī)數(shù)。t時(shí)刻時(shí)粒子i的加速度可表示為

(7)

式中，Mii(t)為t時(shí)刻粒子i的慣性質(zhì)量。粒子i在d維空間中的移動(dòng)速度以及位置的更新可以表示為

Vd，i(t+1)=rand×Vd，i(t)+accd，i(t)

(8)

pd，i(t+1)=pd，i(t)+Vd，i(t+1)

(9)

式中，Vd，i(t)，pd，i(t)分別為第t次粒子移動(dòng)速度和位置。每個(gè)粒子的質(zhì)量是通過第t次的個(gè)體適應(yīng)度計(jì)算得到，可表示為

(10)

式中：fiti(t)為第t時(shí)刻個(gè)體的適應(yīng)度；worst(t)和best(t)分別為第t時(shí)刻種群中個(gè)體最差和最好的適應(yīng)度。假設(shè)個(gè)體的引力質(zhì)量、慣性質(zhì)量和個(gè)體的質(zhì)量相等，則有

Mai=Mpi=Mii=Mi,i=1,2,…,N

(11)

(12)

可知當(dāng)粒子的適應(yīng)度較優(yōu)時(shí)，粒子的質(zhì)量就較大，自身的移動(dòng)速度就會(huì)較慢，會(huì)不斷吸引適應(yīng)度較差的粒子向適應(yīng)度更好的位置移動(dòng)。

3 改進(jìn)GSA

在GSA的過程中，當(dāng)兩個(gè)粒子之間的距離較近時(shí)，粒子就可能陷入局部最優(yōu)解并且很難跳出局部最優(yōu)值，從而降低了搜索效率。本文結(jié)合遺傳算法提出一種改進(jìn)的GSA。首先通過對(duì)初始種群適應(yīng)度計(jì)算，按照適應(yīng)度的大小對(duì)初始種群中的粒子進(jìn)行排序。然后將初始種群分為優(yōu)解組和劣解組兩組，對(duì)優(yōu)解組進(jìn)行隨機(jī)交叉的操作，對(duì)劣解組進(jìn)行變異的操作；將交叉變異得到的新種群與當(dāng)前迭代中的初始種群進(jìn)行重新排序，取適應(yīng)度值最好的前一半種群作為本次迭代后的更新種群；并且引入了概率控制，在整個(gè)算法初期，主要進(jìn)行全局搜索，加大搜索范圍，在后期則進(jìn)行局部挖掘，增加搜索結(jié)果的精度，通過控制進(jìn)行交叉變異操作的概率，能夠有效地發(fā)揮全局搜索和局部挖掘的能力。

首先生成種群數(shù)量為N的初始種群，計(jì)算所有粒子的適應(yīng)度值fitness(t)并且進(jìn)行排序，然后根據(jù)適應(yīng)度將適應(yīng)度較好的一半粒子放入到優(yōu)解組當(dāng)中，再將另一半較差的粒子放入到劣解組當(dāng)中，將優(yōu)解組進(jìn)行遺傳交叉操作，用優(yōu)解組的父代種群產(chǎn)生優(yōu)良的子代種群，并且在本改進(jìn)算法中進(jìn)行交叉操作的兩個(gè)粒子是隨機(jī)選取，可表示為

(13)

劣解組采用變異操作，通過變異操作產(chǎn)生優(yōu)良的子代，充分利用每個(gè)個(gè)體；變異可表示為

(14)

對(duì)改進(jìn)算法流程描述如下，改進(jìn)算法的流程圖，如圖1所示。

圖1 算法流程圖

步驟1設(shè)置算法初始參數(shù)。

步驟2隨機(jī)生成初始種群。

步驟3計(jì)算初始種群的適應(yīng)度，初始化個(gè)體最優(yōu)值和全局最優(yōu)值。

步驟4根據(jù)概率控制判斷是否進(jìn)入遺傳交叉和變異，若rand>p則進(jìn)入步驟5，否則進(jìn)入步驟7。

步驟5進(jìn)行分組策略，根據(jù)適應(yīng)度優(yōu)劣排列粒子順序，將適應(yīng)度較好的一半粒子放入優(yōu)解組進(jìn)行遺傳交叉操作，將適應(yīng)度較差的另一半粒子放入劣解組進(jìn)行變異操作。

步驟6進(jìn)行精英策略，將優(yōu)解組產(chǎn)生的子代和劣解組產(chǎn)生的子代合并為一組，再與當(dāng)前迭代次數(shù)下的初始種群建立一個(gè)擴(kuò)大種群，將擴(kuò)大種群按適應(yīng)度排序，選取適應(yīng)度較好的一半種群作為精英種群參加下一輪迭代。

步驟7用式(8)、式(9)更新粒子的速度與位置。

步驟8更新粒子個(gè)體最優(yōu)值，全局最優(yōu)值。

步驟9判斷是否滿足迭代終止標(biāo)準(zhǔn)，若算法迭代次數(shù)達(dá)到最大迭代次數(shù)T，則終止算法，否則跳回步驟3。

4 算法仿真測(cè)試

4.1 參數(shù)設(shè)置

為了驗(yàn)證改進(jìn)GSA的效果，基于測(cè)試函數(shù)對(duì)粒子群算法、GSA和改進(jìn)GSA進(jìn)行對(duì)比分析。其中粒子群(particle swarm optimization，PSO)算法是由Kennedy等[20]提出的模擬鳥群覓食行為的一種群體協(xié)作優(yōu)化方法。在參數(shù)設(shè)置上，為保證各算法可比性，3個(gè)算法種群規(guī)模均取30，函數(shù)維度均取10，最大迭代次數(shù)取1 000。對(duì)于PSO，其學(xué)習(xí)因子c1和c2分別為1.5和0.5，慣性常數(shù)取0.726，上述參數(shù)是標(biāo)準(zhǔn)PSO的默認(rèn)參數(shù)，已被證明適用于各類優(yōu)化問題；對(duì)于GSA和本文算法，其萬有引力系數(shù)G0皆取10，下降系數(shù)α皆取20。

4.2 測(cè)試函數(shù)

選取5個(gè)函數(shù)標(biāo)準(zhǔn)測(cè)試函數(shù)，將本文改進(jìn)的算法與PSO和GSA算法進(jìn)行比較。其中f1為單峰函數(shù)，主要用來測(cè)試算法的精度與收斂速度；f2為多峰函數(shù)。當(dāng)兩個(gè)函數(shù)維度取2時(shí)，函數(shù)圖像如圖2所示。兩個(gè)函數(shù)用來測(cè)試算法的收斂速度和精度，以及跳出局部最優(yōu)解和全局搜索的能力。兩個(gè)測(cè)試函數(shù)的全局最優(yōu)值均為xi=0,(i=1∶n),f(x)=0。為了消除隨機(jī)性可能對(duì)算法造成的影響，每個(gè)測(cè)試函數(shù)都用3種算法運(yùn)行了100次，試驗(yàn)所處環(huán)境為：MATLAB R2018a，CPU為i5-9400F@2.90 GHz。以均值和方差作為算法的評(píng)判標(biāo)準(zhǔn)。均值表示算法的精度，方差則能反映得到結(jié)果的穩(wěn)定性，收斂曲線能反映算法的收斂速度。f1為單峰值球函數(shù)，所以結(jié)果越接近0就表示越精確，故設(shè)置成功率來顯示其優(yōu)化的效率，當(dāng)適應(yīng)度小于1×10-20，就認(rèn)為搜索成功；f2函數(shù)在最優(yōu)值附近還有許多局部極小值，故當(dāng)適應(yīng)度分別小于1×10-10則認(rèn)為搜索成功。

(a) Sum of different power函數(shù)圖像

(1) Sum of different power函數(shù)：x∈[-1,1]

(15)

(2) Ackley’s Path函數(shù)：x∈[-1.5,1.5]

(16)

4.3 仿真結(jié)果分析

為了比較改進(jìn)GSA，PSO和GSA算法的優(yōu)劣，圖3列出了兩種測(cè)試函數(shù)的平均適應(yīng)度曲線圖，各算法100次搜索的統(tǒng)計(jì)結(jié)果，如圖3所示。從平均適應(yīng)度曲線的數(shù)據(jù)可以看出，改進(jìn)的優(yōu)化算法在精度上都要優(yōu)于PSO和GSA算法，改進(jìn)算法收斂時(shí)的結(jié)果對(duì)于f1測(cè)試函數(shù)都要小于1×10-20，對(duì)于測(cè)試函數(shù)f2其收斂時(shí)的精度也明顯高于另外兩種算法。

(a) f1

對(duì)于f1函數(shù)改進(jìn)的優(yōu)化算法得到的適應(yīng)度相比其他算法其探索點(diǎn)的精度顯著提高，在優(yōu)化的成功率上，其達(dá)到數(shù)量級(jí)1×10-15以下的成功率為100%遠(yuǎn)遠(yuǎn)高于其他優(yōu)化算法。對(duì)于測(cè)試函數(shù)f2改進(jìn)算法探索到最優(yōu)位置的成功率接近100%。且測(cè)試函數(shù)f2其全局最優(yōu)值相對(duì)于局部極小值的適應(yīng)度相差很小且局部最優(yōu)值數(shù)量較多，因此很難找到最優(yōu)位置，改進(jìn)GSA的成功率有98%，其他算法則很難找到最優(yōu)值；這表明改進(jìn)GSA顯著提升了跳出局部最優(yōu)解的能力。將3種算法相對(duì)于兩種測(cè)試函數(shù)的計(jì)算效率進(jìn)行了對(duì)比，列出了完成一次尋優(yōu)操作所需時(shí)間具體時(shí)間，如表1所示。可以看到，相比基本GSA和PSO算法，由于改進(jìn)GSA算法加入了交叉變異等操作；因此，其效率上不占優(yōu)勢(shì)，但也未顯著增加計(jì)算時(shí)間。

表1 計(jì)算效率對(duì)比

5 南中環(huán)橋模型修正

5.1 工程概況

南中環(huán)橋是連接太原市東西城區(qū)的重要橋梁，如圖4所示。該橋主橋結(jié)構(gòu)體系為組合式系桿拱、鋼管混凝土蝶形拱橋，跨徑布置為(60+180+60) m。主橋跨中為鋼混疊合梁段，全長(zhǎng)92 m，主梁高為3 m；主橋邊跨60 m及中跨44 m梁段，采用的是混凝土箱梁，梁高為3.0～6.5 m。主橋標(biāo)準(zhǔn)斷面寬度53 m，跨中寬度為71 m。主跨是由4根主副拱肋組合而成的下承式鋼管混凝土拱橋，主拱肋外傾16°，矢跨比為1/4.326；副拱肋外傾26.882°；主拱肋和副拱肋之間利用拉桿連接。

圖4 南中環(huán)橋

5.2 荷載試驗(yàn)概況

南中環(huán)橋荷載試驗(yàn)分為靜載試驗(yàn)和動(dòng)載試驗(yàn)兩個(gè)部分，靜載試驗(yàn)主要測(cè)試橋跨結(jié)構(gòu)在試驗(yàn)車輛靜載作用下的位移；動(dòng)載試驗(yàn)主要測(cè)試結(jié)構(gòu)的固有振動(dòng)特性和沖擊系數(shù)。限于篇幅，僅介紹與模型修正相關(guān)的試驗(yàn)工況。模型修正共選擇了兩個(gè)靜載試驗(yàn)工況(A-A、B-B)，測(cè)試在車輛下主梁和主拱關(guān)鍵截面的位移。圖5為橋梁立面圖及位移測(cè)點(diǎn)布置圖。靜載試驗(yàn)車輛為40 t三軸載重汽車，軸距為(1.8+3.3) m。動(dòng)力測(cè)試部分，為獲取橋梁的固有振動(dòng)特性，利用加速度傳感器測(cè)試主梁在環(huán)境脈動(dòng)荷載下的三向加速度響應(yīng)。在主梁縱向平均每14 m布置一個(gè)加速度測(cè)點(diǎn)，采樣頻率為80 Hz。利用測(cè)試得到的加速度響應(yīng)，結(jié)合模態(tài)識(shí)別算法可以獲得橋梁的頻率、阻尼比和振型。靜力位移、頻率的試驗(yàn)值和計(jì)算值的對(duì)比結(jié)果，如表2和表3所示。從表2和表3可知，位移最大相對(duì)誤差達(dá)28.19%，頻率的最大相對(duì)誤差達(dá)到9.51%；其中，計(jì)算值是根據(jù)Midas Civil初始有限元模型計(jì)算得到的。該模型按設(shè)計(jì)圖紙建立，全橋共有8 322個(gè)單元，6 118個(gè)節(jié)點(diǎn)，其中主副吊桿和斜拉桿采用桁架單元，其余部分皆用梁?jiǎn)卧M。主梁采用梁格法模擬，鋼管混凝土拱采用聯(lián)合截面模擬；吊桿張拉時(shí)考慮幾何剛度貢獻(xiàn)，邊界條件按設(shè)計(jì)圖紙?jiān)O(shè)置。

表2 兩種工況下各測(cè)試截面的位移計(jì)算值和試驗(yàn)值

表3 橋面系實(shí)測(cè)自振頻率及其計(jì)算頻率

圖5 橋梁立面圖及位移測(cè)點(diǎn)布置圖(cm)

主梁前兩階豎向振型計(jì)算值和試驗(yàn)值的對(duì)比，如圖6所示。對(duì)應(yīng)的模態(tài)置信值(modal assurance criterion, MAC)分別為0.946和0.963；這表明試驗(yàn)振型與計(jì)算振型吻合良好，反映試驗(yàn)數(shù)據(jù)可靠。上述分析表明：初始有限元模型對(duì)結(jié)構(gòu)真實(shí)剛度分布模擬不夠準(zhǔn)確，需要通過模型修正，從而更好地代表實(shí)際結(jié)構(gòu)行為。

(a) 主梁二階振型

5.3 Kriging模型建立

為了提高模型修正迭代效率，利用Kriging模型替代有限元模型預(yù)測(cè)計(jì)算值。首先選擇修正參數(shù)，然后利用拉丁超立方抽樣獲得參數(shù)在設(shè)計(jì)空間的分布，最后建立Kriging模型并評(píng)估其精度。

在修正參數(shù)選擇方面，因?yàn)樵摌驗(yàn)轭A(yù)制構(gòu)件所以截面誤差一般較小，故截面慣性矩不作為修正參數(shù)。在建立有限元模型時(shí)對(duì)于整個(gè)橋梁進(jìn)行了簡(jiǎn)化，并沒有建立節(jié)點(diǎn)板、螺栓等構(gòu)件，所以在結(jié)構(gòu)的質(zhì)量分布上會(huì)存在影響，而一些細(xì)部的加勁肋則會(huì)影響到模型的整體剛度，因此將參數(shù)的選擇范圍定在各個(gè)構(gòu)件的質(zhì)量密度和彈性模量上。本文分別選擇了鋼管拱主拱、鋼管拱副拱、鋼混疊合梁段、混凝土段和挑臂5個(gè)構(gòu)件的彈性模量和質(zhì)量密度共10個(gè)參數(shù)進(jìn)行敏感性分析。結(jié)合敏感性分析結(jié)果，最終選取鋼管拱主拱彈性模量E1、鋼管拱主拱密度ρ1、鋼混疊合梁段彈性模量E2、鋼混疊合梁段密度ρ2、混凝土段主梁彈性模量E3和混凝土段主梁質(zhì)量ρ3這6個(gè)參數(shù)作為待修正參數(shù)，其初始值分別為2.06×108kN/m2，78.5 kN/m3，2.06×108kN/m2，78.5 kN/m3，3.45×107kN/m2，26 kN/m3。

采用拉丁超立方抽樣法對(duì)選定的修正參數(shù)進(jìn)行試驗(yàn)設(shè)計(jì)，對(duì)6個(gè)設(shè)計(jì)參數(shù)進(jìn)行了81次抽樣，參數(shù)的抽樣范圍設(shè)定在±30%以內(nèi)，再將設(shè)計(jì)得到的修正參數(shù)分別代入初始有限元模型中得到對(duì)應(yīng)工況下的位移以及頻率響應(yīng)。以樣本點(diǎn)作為輸入，節(jié)點(diǎn)位移值及模態(tài)頻率響應(yīng)值作為輸出構(gòu)建Kriging模型。選取鋼管主拱彈模E1、鋼混疊合梁段彈模E2兩個(gè)設(shè)計(jì)參數(shù)和該橋有限元模型一階頻率繪制三維Kriging響應(yīng)面圖。一階頻率對(duì)于修正參數(shù)E1和E2的Kriging響應(yīng)面模型以及對(duì)應(yīng)的均方誤差(mean square error，MSE)，如圖7所示。其中曲面為Kriging模型響應(yīng)面，散點(diǎn)為樣本點(diǎn)，可知樣本點(diǎn)基本位于曲面上。從對(duì)應(yīng)的MSE圖可以看出，誤差基本都不超過3×10-9，說明Kriging精度較高，可替代有限元模型進(jìn)行結(jié)構(gòu)響應(yīng)預(yù)測(cè)。

圖7 一階頻率與E1、E2之間函數(shù)擬合圖

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

考慮到實(shí)際結(jié)構(gòu)損傷位置是未知的，為了證明模型修正方法的有效性，首先利用南中環(huán)橋有限元模型進(jìn)行數(shù)值驗(yàn)證。給參數(shù)E1,E2,E3設(shè)定數(shù)值損傷，其中：相較于初始值，參數(shù)E1,E2設(shè)定+10%的數(shù)值損傷，參數(shù)E3設(shè)定-10%的數(shù)值損傷，以未設(shè)損傷的模型代表初始有限元模型，施加損傷后的模型為真實(shí)的有限元模型。以3個(gè)設(shè)計(jì)參數(shù)作為輸入，以前6階頻率作為輸出，建立Kriging模型，然后以本文中的3種優(yōu)化算法對(duì)建立的Kriging模型進(jìn)行尋優(yōu)，PSO算法仍取其默認(rèn)參數(shù)；GSA和改進(jìn)GSA算法的萬有引力系數(shù)取為10，下降系數(shù)為20；各算法的種群規(guī)模均為30，最大迭代次數(shù)為1 000。其設(shè)計(jì)參數(shù)和頻率響應(yīng)的結(jié)果,如表4、表5所示。由表可知，改進(jìn)算法3個(gè)設(shè)計(jì)參數(shù)的修正偏差修正前最高為11.11%，修正后誤差都降到約0.1%；頻率的響應(yīng)值在修正前最大偏差達(dá)到1.78%，修正之后改進(jìn)算法的最大誤差僅為-0.024%。采用本文的改進(jìn)算法結(jié)合Kriging模型的模型修正方法是有效的。其中改進(jìn)GSA、GSA和PSO算法的計(jì)算時(shí)間分別為0.468 s、0.353 s和0.146 s。對(duì)于橋梁模型修正而言，每次迭代中最為耗時(shí)的部分是評(píng)估目標(biāo)函數(shù)，即計(jì)算設(shè)計(jì)參數(shù)修改后的目標(biāo)函數(shù)值，而算法所耗時(shí)間相對(duì)較少；因此，改進(jìn)GSA的計(jì)算效率可以接受。

表4 參數(shù)修正前后對(duì)比

表5 頻率修正前后對(duì)比

5.5 基于實(shí)測(cè)數(shù)據(jù)的模型修正

以位移殘差和前4階頻率殘差構(gòu)建目標(biāo)函數(shù)，分別利用改進(jìn)算法、GSA和PSO對(duì)目標(biāo)函數(shù)進(jìn)行尋優(yōu)，算法在參數(shù)設(shè)置上，PSO算法仍取其默認(rèn)參數(shù)；GSA和改進(jìn)GSA算法的萬有引力系數(shù)取為10，下降系數(shù)為20；各算法的種群規(guī)模均為30，最大迭代次數(shù)為1 000。修正后的參數(shù)、頻率和位移結(jié)果分別如表6、表7和表8所示。

由表6可知，參數(shù)鋼管拱主拱彈性模量E1、鋼混疊合梁段彈性模量E2、混凝土段彈性模量E3利用改進(jìn)萬有引力搜索算法的修正結(jié)果相較于初始值分別提高了13.27%，21.72%，24.37%，主要原因?yàn)樵诮⒊跏寄Ｐ蜁r(shí)對(duì)橋梁整體進(jìn)行了簡(jiǎn)化，在建立有限元模型時(shí)并沒有考慮局部加勁肋和橫隔板，并且鋼管拱以及鋼混疊合梁段在實(shí)際橋梁中數(shù)量較多，所以有限元模型比實(shí)橋的剛度較低；鋼管拱主拱密度ρ1的修正值都與實(shí)際值差別不大，都略微偏小；而鋼混疊合梁段密度ρ2的修正值都比實(shí)際值偏大，混凝土段質(zhì)量ρ3的修正值則比實(shí)際值都較為偏小，原因是在混凝土施工質(zhì)量存在一定的離散性，使分布不均勻，并在模型修正中體現(xiàn)出來。需要指出的是，實(shí)際橋梁參數(shù)修正前后的變化參數(shù)修正前后的變化，僅表示模型整體剛度和質(zhì)量分布與實(shí)際情況的差異，并非材料參數(shù)發(fā)生了變化，即修正后的參數(shù)可理解為等效彈性模量和等效質(zhì)量密度。并且通過查看設(shè)計(jì)圖紙可以得知鋼混疊合梁段鋼梁節(jié)段一的總質(zhì)量為107 434.7 kg，而通過Midas Civil查詢單元質(zhì)量可知初始有限元模型鋼梁節(jié)段一的總質(zhì)量為103 025.2 kg，即有限元模型中的質(zhì)量低估了約4.28%，這與鋼混疊合梁段密度ρ2的修正值3.77%較為接近。

表6 參數(shù)修正前后對(duì)比

由表7可知，利用改進(jìn)GSA、GSA和PSO修正后頻率值的對(duì)比，其相對(duì)誤差皆為修正前后計(jì)算值相對(duì)于實(shí)測(cè)值的誤差。從整體看南中環(huán)橋的實(shí)測(cè)頻率比計(jì)算頻率都偏大，說明有限元模型的整體剛度都偏低，與參數(shù)的修正結(jié)果符合，修正后的頻率值都有所提升，說明修正結(jié)果可作參考。相比初始模型，各算法修正后的頻率相對(duì)誤差均有不同程度降低，改進(jìn)GSA取得了相對(duì)于PSO和GSA更好的修正結(jié)果。

表7 頻率修正前后對(duì)比

由表8可知，3個(gè)優(yōu)化算法修正后位移值的對(duì)比，除了節(jié)點(diǎn)1和節(jié)點(diǎn)3在修正前其計(jì)算值和位移值差距不大導(dǎo)致修正空間較小，其余7個(gè)節(jié)點(diǎn)都相對(duì)于修正前精度有明顯的提高。整體上有限元模型計(jì)算值與實(shí)測(cè)值相比都偏大，也證明了有限元模型在整體剛度上都偏小，修正后結(jié)果都比修正前更加接近實(shí)測(cè)值。修正前，實(shí)測(cè)值與計(jì)算值的最大誤差達(dá)到了28.19%，在修正后，改進(jìn)算法所得結(jié)果的相對(duì)誤差為13.62%，在精度上得到明顯的提升，且GSA和PSO的修正結(jié)果相對(duì)誤差分別為16.39%和15.90%，改進(jìn)算法也比這兩種算法修正精度要高。上述結(jié)果表明：經(jīng)過模型修正，模型的位移值計(jì)算值和測(cè)試值更為接近，其相對(duì)誤差減小，表明了修正后更能反映實(shí)際情況。相對(duì)于GSA和PSO，改進(jìn)GSA的修正精度更高。

表8 位移修正前后對(duì)比

6 結(jié) 論

通過試驗(yàn)與有限元模擬，基于Kriging模型對(duì)南中環(huán)橋進(jìn)行了有限元模型修正。為提高修正的精度，提出了一種改進(jìn)的優(yōu)化算法，得出了以下結(jié)論：

(1) 對(duì)于測(cè)試函數(shù)的修正結(jié)果，改進(jìn)的優(yōu)化算法對(duì)于普通GSA和PSO算法都有著較高的成功率和穩(wěn)定性，從而提高了優(yōu)化結(jié)果合理解的可能性；對(duì)于南中環(huán)橋的修正結(jié)果，改進(jìn)的優(yōu)化算法獲得更小的目標(biāo)函數(shù)，且修正后頻率和位移整體相對(duì)誤差較小，剛度分布更為合理。

(2) 經(jīng)過修正后位移和頻率的計(jì)算值與試驗(yàn)值之間的誤差更小，其位移頻率除個(gè)別值外其誤差都有顯著降低，其中位移較大的兩個(gè)節(jié)點(diǎn)誤差分別從23.92%和28.19%降低到14.35%和13.62%，表明修正后的模型能更好代表實(shí)際結(jié)構(gòu)。

(3) 因?yàn)閷?shí)際結(jié)構(gòu)的復(fù)雜性和試驗(yàn)結(jié)果的不確定性，故模型修正的目標(biāo)函數(shù)在修正參數(shù)的空間分布較復(fù)雜，本文提出的優(yōu)化算法對(duì)修正結(jié)果的精度和穩(wěn)定性有所提升。