基于全局敏感性的自適應模型更新混合模擬方法

陳再現, 王紀偉,2, 鐘煒彭, 劉 鋮, 李明剛

(1.哈爾濱工業大學(威海) 海洋工程學院, 山東 威海 264209; 2.中建三局第一建設工程有限責任公司,湖北 武漢 430040)

混合模擬把計算機難以精確模擬的復雜非線性部分作為物理子結構進行真實的物理加載,其余部分作為數值子結構在計算機中建模計算,在20世紀70年代初由日本學者Hakuno等[1]提出,并在20世紀70年代中期由日本學者Takanashi等[2]第一次實現。混合模擬減少了試驗成本和試驗場地及設備的要求,并且也較為真實地反映了非線性較強構件的受力情況。因此被普遍認為是一種經濟高效且極具有發展前景的新型抗震試驗技術。

數值子結構的精確程度是影響混合模擬結果的重要因素。模型更新依據試驗子結構反饋數據不斷修正數值子結構,因此能較好地提升數值子結構的模型精度。近年來,模型更新備受眾多學者關注。Yang等[3]首次將模型更新引入到混合模擬,利用神經網絡對數值結構建模,通過試驗子結構反饋數據不斷修正神經網絡,并在二層剪切模型中通過數值模擬和真實試驗進行了驗證。Elanwar等[4-5]以數值子結構和物理子結構的應力差值作為目標函數,對材料本構參數進行識別,數值模擬結果表明,基于材料本構模型參數的更新方法能提高混合模擬的精度。Yang等[6-7]分別開展了鋼框架和鋼筋混凝土橋墩的混合模擬工作,從材料本構模型層面對參數進行更新并驗證了其對混合模擬精度的提升。Zhong等[8]提出采用均勻設計的方法對鋼結構和鋼筋混凝土結構材料參數進行識別,該方法極大地降低了樣本數量,提高了識別效率,該方法在二層框架中以數值模擬的方式得到了驗證。王濤等[9-10]提出統計容積卡爾曼濾波器的模型更新方法,并引入多尺度模型更新概念,通過二層帶有自復位摩擦耗能支撐框架結構的數值模擬分析驗證了方法的可行性。

模型更新首先依據尋優方法創造樣本空間,然后利用試驗子結構反饋數據確定最優樣本點,最后更新至數值結構。從模型更新的過程來看,它是求解逆問題,根據試驗子結構反饋數據識別多參數,樣本空間的選取不當直接影響參數識別的效果。陳再現等[11-12]提出基于全局敏感性(global sensitivity, GS)的混合模擬試驗方法,通過參數敏感度系數確定參數的搜索范圍即樣本空間能使參數更容易收斂到有效值。然而,在多參數識別過程中,僅通過固定的靈敏度系數去指導模型更新易導致參數陷入局部最優,從而降低模型更新的精度。基于此,將動量和自適應步長的概念引入模型更新中,提出基于全局敏感性的自適應(global sensitivity adaptive, GS-A)混合模擬試驗方法在一定程度上彌補了僅用敏感度系數確定樣本空間的缺陷,提升了參數收斂的速度和模型更新的精度。下文通過介紹所提出方法的基本原理、用以驗證該方法的數值模型和更新參數的選擇以及數值模擬方案,通過33個數值模擬算例驗證本文方法的有效性。

1 GS-A混合模擬模型更新方法基本原理及流程

1.1 基本原理

陳再現等推出的GS模型更新混合模擬方法使用均勻設計選擇樣本空間時,以上一步最優參數為搜索中心,各參數搜索步長與對應的敏感性系數成正比,利用相應的均勻設計表[13]生成驗算子結構樣本空間。

然而在參數識別過程中,若固定中心點(上一步最優參數組合),創建樣本空間,對當前步最優參數進行搜索,會出現與參考模型差異較大的參數不能迅速被識別,設置落入局部最優點,因此基于全局敏感性的自適應(GS-A)模型更新混合模擬方法引入動量概念,即參考歷史參數變化,引導參數朝著最優值更快收斂,動量表達式為

(1)

在參數識別的最初階段,參數組合距離目標函數的最優解距離很遠,但是隨著參數識別的次數變大,參數組合會越來越接近最優解,因此本文提出對搜索步長進行自適應處理,依據每個參數的變化量大小來調整各個參數的步長,使參數組合在接近最優解時,步長將適當縮小,使參數組合更加趨近于最優解。自適應步長計算公式見式(2)～式(3)

(2)

(3)

式中,ε為維持數值穩定性所添加的常數,本文取為10-8。

(4)

(5)

式中,TSi為t時刻,第i個參數的總階靈敏度系數。

1.2 基本流程

數值模擬所采用的數值模型[14]由全結構數值模型、物理子結構和驗算子結構三部分組成,如圖1所示上部三個部分。將地震激勵輸入到全結構數值模型,通過時程分析獲取當前步的內力和位移,并輸入到物理子結構和驗算子結構中,獲得恢復力,不同的是,物理子結構的恢復力為參考目標,驗算子結構是物理子結構的數值模擬空間(通過優化算法選擇),兩者除更新參數外均相同,依據目標函數選取與物理子結構恢復力誤差最小的驗算子結構,將其參數更新到計算全結構中,從而進行下一地震步的計算。如圖1所示模型更新基本流程,具體步驟如下:

圖1 模型更新流程圖

步驟1利用動力時程分析求出數值全結構結構響應,提取出子結構劃分界面處的位移及內力(界面位移yps,軸力Nps和彎矩Mps),并發送到物理子結構和驗算子結構。

步驟2物理子結構和驗算子結構接收到位移與內力,完成加載后得到各自的恢復力,并根據誤差判斷規則,選擇是否啟動更新程序,目標函數選取與物理子結構恢復力誤差最小的驗算子結構,如式(6)所示,如果啟動,進入步驟3,如果不啟動,返回步驟1進行下一步動力時程分析。

步驟3依據式(4)所確定的參數搜索范圍,利用均勻設計得到驗算子結構樣本空間,識別程序比較樣本空間中多個驗算子結構與物理子結構恢復力的相對累積誤差,選出誤差最小的驗算子結構,將其參數組合作為更優參數組合發送到全結構數值模型進行參數更新,并進行下一步動力時程分析。

步驟4直到地震動記錄結束,終止程序,否則,重復步驟1～步驟4。

(6)

式中:Rps(i)為物理子結構恢復力,作為參考值;Rys(i)為驗算子結構恢復力。

2 數值模擬方案

2.1 基本模型構建

選用一榀二層一跨的鋼框架框結構,一層和二層高度分別為3 600 mm和3 000 mm,跨度6 000 mm,框架幾何尺寸、質量及荷載分布如圖2所示,質量及荷載分布參照GB 50011—2010《建筑抗震設計規范》與GB 50009—2012《建筑結構荷載規范》計算所得,截面信息分別如圖2和圖3所示。

圖2 框架模型

圖3 鋼框架截面

物理子結構與數值子結構的界面設置在彈性反彎點處[15],由D值法計算得到鋼框架的反彎點分別為底柱長度的0.58倍處,如圖2所示。采用OpenSees有限元軟件進行數值建模,鋼材本構模型采用雙線性模型(Steel01),模型假定基礎與地基剛性連接,框架結構采用基于力的非線性梁柱單元(nonlinear beam-column element),每個單元選取5個Gauss-Lobatto,截面為纖維截面。

2.2 更新參數確定及模擬方案

選用雙線性鋼材本構模型(Steel01),更新參數為8個:彈性模量Es,屈服強度fys以及硬化系數bs,尺寸參數為高度h、寬度d、腹板厚度t1和翼緣厚度t2,以及柱長ls。參數的參考值設置如表1所示,將各個參數比例化,設定參數的參考值為1,模型的材料誤差參數浮動范圍設為(0.7～1.3),尺寸誤差參數浮動范圍設為(0.98～1.02)。

表1 鋼結構參數參考值

表2 初始參數組

為驗證GS-A方法模型更新效果,把不更新混合模擬(hybrid simulation,HS)和GS模型更新方法作為對照組。HS用于反映傳統混合模擬中各參數是帶有誤差的統計值的情況并且在試驗過程中參數不被更新。GS方法和GS-A方法分別為在混合模擬中采用基于全局敏感性的模型更新方法和基于全局敏感性的自適應模型更新方法,在試驗過程中參數不斷地被更新。數值模擬方案如表3所示。

表3 數值模擬方案

3 數值模擬結果

各參數對一層位移相對累計誤差的敏感性在地震的各個時間段是不斷變化的,因此對試驗模型參數進行分階段敏感性分析,并以此結果作為模型更新驗算子結構樣本空間選擇的依據,敏感性分析方法采用Sobol法[16]。設置地震峰值加速度( peak ground acceleration,PGA)為510 gal的El Centro波作為地震激勵,取持續6 s的地面運動記錄,積分步長取0.01 s。地震分階段靈敏度系數變化如圖4所示。根據Sobol分析結果,在結構處于彈性階段時(2 s之前),參數彈性模量和翼緣寬度的敏感度系數最大,在結構經歷塑性階段時(2 s之后),參數屈服強度敏感度系數開始增大,屈服強度、彈性模量和翼緣厚度的敏感性排名前三。

圖4 地震分階段鋼框架模型參數靈敏度系數變化圖

參考解(Ref)是模型參數為參考值時的時程分析結果。對照組分別為參數不更新時混合模擬(HS)的結果和采用基于全局敏感性(GS)的模型更新方法進行模型更新的結果。采用物理子結構位移相對累積誤差和峰值位移相對誤差作為模擬結果精度的評價指標,分別如式(7)和式(8)。

(7)

式中:yref(i)為物理子結構位移參考值,由時程分析得出;yps(i)為物理子結構位移。

(8)

式中:y_maxref為物理子結構峰值位移參考值,由時程分析得出;y_maxps為物理子結構峰值位移。

鋼框架模型更新效果如圖5所示,在11個初始誤差參數組中,相對于HS,使用GS-A方法和GS方法更新完成后模型的物理子結構處位移時程與參考位移時程的相對累積誤差均有所降低,并且,兩種模型更新方法對峰值位移相對誤差均有所改善,但是GS-A更新效果明顯優于GS,GS-A物理子結構處位移的相對累積和峰值位移相對誤差均在5.7%和2.7%以內,而GS在案例7中甚至達到了52.9%和28.4%;11個初始誤差參數組中,GS-A方法和GS方法的物理子結構處位移的相對累積誤差和峰值位移相對誤差的均值和提升效果如表4所示,GS-A方法物理子結構處位移的相對累積和峰值位移相對誤差的均值為2.3%和1.3%,均優于GS方法的17.8%和9.4%,相對于HS,分別提升了97.1%和96.6%,提升效果均優于GS方法的78.3%和74.5%。

表4 模型更新提升效果

圖5 鋼框架(steel01本構模型)更新效果圖

GS方法和GS-A方法對11個初始誤差參數組的參數識別結果的平均誤差如表5所示,表中按參數敏感性由大到小排序。可見,GS-A方法的參數識別誤差相對于GS方法有所降低,值得注意的是,對于敏感性最小的硬化系數bs,GS-A方法能有效識別,且平均誤差僅有1.3%,GS方法則識別失效,平均誤差達到27.9%。

表5 模型更新后參數平均誤差

GS方法和GS-A方法在模型更新后段(4～6 s)中各參數的平均相對誤差和標準差如表6所示,表中按參數敏感性由大到小排序。GS-A方法參數平均相對誤差均小于2.3%,而GS方法對敏感性較大的參數彈性模量的識別平均誤差達到4.9%,對于敏感性較小的硬化系數平均誤差達到27.9%。對于參數識別過程的標準差,GS-A方法均小于GS方法,因此GS-A方法的穩定性明顯優于GS方法。

表6 模型更新(4～6 s)參數平均誤差和標準差

在11個案例中,案例9所對應的GS方法的位移的相對累積誤差和峰值相對誤差均比HS大,以此選取案例9的參數識別過程、位移時程與滯回曲線對比示例。

案例9的參數識別過程對比如圖6所示,GS方法和GS-A方法在模型更新后段(4～6 s)中各參數的標準差如表7所示,GS-A方法的標準差均低于GS方法,因此GS-A方法的參數識別的穩定性優于GS方法。值得注意的是,在GS方法中,各參數在識別過程中不斷上下波動,說明參數陷入局部最優值,從而來回震蕩,反觀GS-A方法,動量使樣本空間中心向參數變化方向移動,參數更快地收斂,自適應系數調整樣本空間的范圍,防止參數組合陷入局部最優。GS方法和GS-A方法參數最終識別值與參考值的比值如表8所示。由表8可知,GS-A方法的參數識別誤差相對于GS方法均有所降低,提升了模型更新精度。

表7 模型更新(4～6 s)參數標準差(案例9)

表8 參數最終識別值與參考值比值(案例9)

圖6 參數更新效果圖(案例9)

案例9的物理子結構位移時程和滯回曲線對比見圖7和圖8,GS-A方法的位移時程曲線和滯回曲線相對HS和GS方法更貼近于參考解(Ref)。GS-A方法的物理子結構位移與參考解的相對累積誤差為2.27%,遠小于HS的23.05%和GS方法的47.16%,因此GS-A方法顯著地提升了模型更新效果。

圖7 位移曲線對比圖(案例9)

圖8 滯回曲線對比圖(案例9)

4 結 論

(1)提出了GS-A模型更新混合模擬方法,在鋼框架中進行了11種不同初值參數誤差組的數值模擬用以驗證該方法的適用性,驗證了動量和自適應策略的可行性。模擬結果顯示,GS-A方法的相對累積誤差和物理子結構峰值位移相對誤差的均值分別為2.3%和1.3%,相對于GS方法的17.8%和9.4%有顯著提升,改善了模型更新效果。

(2)在參數識別精度方面,在GS方法中敏感性較大的參數彈性模量的平均相對誤差達到了6.0%,敏感性較小的參數硬化系數識別效果欠佳,平均相對誤差達到了27.1%,而GS-A方法各參數的平均相對誤差均在1.4%以內,GS-A方法相對于GS方法有效地提升了參數識別的精度。

(3)在參數識別后段(4～6 s),GS-A方法中的各參數標準差平均值均小于GS方法,因此GS-A方法提升了參數識別的穩定性。

(4)案例9中,GS方法參數組合陷入局部最優解,因此導致模型更新效果比不更新還差。GS-A方法解決了參數陷入局部最優解的問題,參數識別效果的精度和穩定性均優于GS方法,且位移時程和滯回曲線相對于HS和GS方法均更貼近于參考解,物理子結構位移相對累積誤差僅有2.3%,顯著地提升了模型更新效果。