基于POD -RBF代理模型的迭代更新反演方法

呂小龍， 黃 丹，姜冬菊

(河海大學(xué) 力學(xué)與材料學(xué)院，南京 211100)

1 引 言

水工大壩等大體積混凝土工程結(jié)構(gòu)的材料參數(shù)對(duì)其健康和安全狀態(tài)的監(jiān)測起著十分重要的作用，但結(jié)構(gòu)內(nèi)部的材料參數(shù)往往難以直接得到，于是通過位移等可觀測量來反演壩體材料參數(shù)便成為一種廣泛使用的途徑。然而，現(xiàn)有的大量大體積結(jié)構(gòu)材料參數(shù)反演方法往往需要通過有限元等數(shù)值方法進(jìn)行大規(guī)模正分析，且正分析耗時(shí)量隨著待反演參數(shù)的增加呈幾何式增長。對(duì)于復(fù)雜的大體積結(jié)構(gòu)，甚至可能導(dǎo)致較多參數(shù)反演問題的正分析環(huán)節(jié)難以施行。通過建立代理模型以取代大量的有限元計(jì)算，可以花費(fèi)很少的時(shí)間獲得足夠精度的響應(yīng)值，因而廣泛應(yīng)用于結(jié)構(gòu)優(yōu)化、可行性分析和預(yù)測反演等問題中[1]。

Bhosekar等[2]對(duì)各種常見代理模型的構(gòu)建、應(yīng)用及優(yōu)勢進(jìn)行了總結(jié)。代理模型本質(zhì)上是輸入和輸出之間的一種映射關(guān)系，常見的映射模型包括多項(xiàng)式模型[3]、Kriging插值模型[4，5]、徑向基函數(shù)模型[6-8]以及人工神經(jīng)網(wǎng)絡(luò)模型[9,10]等。其中采用本征正交分解POD(Proper Orthogonal Decomposition)和徑向基函數(shù)RBF(Radial Basis Functions)插值相結(jié)合構(gòu)建的代理模型，能同時(shí)保證較高的擬合精度和插值精度[11]，適用于大體積結(jié)構(gòu)的模擬預(yù)測。本文研究思路即來自于此。此外，為了在保證反演精度的前提下盡可能降低有限元正分析計(jì)算量，本文結(jié)合粒子群算法的全局尋優(yōu)能力和高斯-牛頓法的局部快速收斂能力，構(gòu)建了一種組合式的迭代更新反演模型，并通過典型大體積混凝土大壩的分區(qū)彈性模量反演驗(yàn)證其反演精度和計(jì)算效率。

2 代理模型構(gòu)造及原理

2.1 本征正交分解(POD)

POD是統(tǒng)計(jì)學(xué)上用來分析數(shù)據(jù)的常用方法[12]，在數(shù)據(jù)科學(xué)領(lǐng)域也稱為主成分分析。以典型水工大壩為例，考慮將大壩分成p個(gè)彈模分區(qū)，設(shè)置n組彈模樣本得到彈模參數(shù)樣本矩陣X=[x1,x2,…,xn]，其中xi(i=1,2,…,n)為包含一組彈模值的列向量。在相同荷載作用下，其分區(qū)彈模不同的取值對(duì)應(yīng)不同的位移響應(yīng)，設(shè)位移響應(yīng)矩陣為Y=[y1,y2,…,yn]，其中yi(i=1,2,…,n)為m維位移響應(yīng)列向量，表示在一組分區(qū)彈模值下大壩上m個(gè)測點(diǎn)的位移測值所組成的列向量，稱之為快照(Snapshot)，n為樣本數(shù)量。POD的關(guān)鍵目的就是找出這些快照的內(nèi)在聯(lián)系，提取數(shù)據(jù)的主要特征，進(jìn)而對(duì)數(shù)據(jù)進(jìn)行降維。對(duì)Y進(jìn)行奇異值分解

Y=Φ·Λ·VT

(1)

式中Φ=[φ1,φ2,…,φm] (φi?Rm,i=1,2,…,m)為Y的左奇異向量,V=[v1,v2,…,vn] (vi?Rn,i=1,2,…,n)為Y的右奇異向量,Λ為m行n列的對(duì)角矩陣，其對(duì)角線上的每個(gè)元素λi(i=1,2,…,min(m,n))為Y的奇異值，即協(xié)方差矩陣YTY的特征值的平方根。對(duì)大部分情況，奇異值λi(i=1,2,…,min(m,n))從大到小衰減的速度很快，前幾個(gè)奇異值的總和通常就能占到所有奇異值總和的90%以上，即

(2)

根據(jù)奇異值分解的性質(zhì)，原數(shù)據(jù)Y在左奇異向量φi(i=1,2,…,n)上的投影是最大的，因此，以前k階左奇異向量為基向量便能夠包含Y的絕大部分信息。由于k通常遠(yuǎn)小于n，故原數(shù)據(jù)在保留絕大部分信息的前提下得到了大幅降維。以前k階左奇異向量為基向量，任意彈模參數(shù)樣本x對(duì)應(yīng)的位移響應(yīng)向量y(x)可以近似表示其線性組合為

(3)

2.2 徑向基函數(shù)(RBF)模型

RBF是一種常用的多變量函數(shù)插值基函數(shù)，其一般形式可記為

(4)

Thin-platesplines函數(shù):

f(r)=r2lnr

(5)

Inverse multiquadric函數(shù)：

f(r)=(r2+c2)-1/2

(6)

以及高斯函數(shù)：

f(r)=exp(-r2/c2)

(7)

式(6,7)的c為一個(gè)大于零的常數(shù)，稱為平滑系數(shù)，其值越大，函數(shù)值變化越平緩，反之越陡峭。綜合已有工作比較，本文采用inverse multiquadric函數(shù)作為徑向基函數(shù)，平滑系數(shù)c=0.5。

將所有已知樣本參數(shù)點(diǎn)和對(duì)應(yīng)的位移響應(yīng)代入插值式(4)中，并寫成矩陣形式

(8)

2.3 形成POD -RBF代理模型

根據(jù)RBF插值方法，可對(duì)方程(3)的組合系數(shù)α(x)進(jìn)行插值，即

(9)

或矩陣形式

α(x)=W·f(r)

(10)

式中W為wi j組成的系數(shù)矩陣，f(r)為f(ri)組成的列向量。將式(10)代入式(3)，并應(yīng)用已有的所有樣本點(diǎn)得到

(11)

(12)

即可求得系數(shù)矩陣W。

由式(3,10,12)可得POD -RBF代理模型

(13)

3 迭代更新反演算法

參數(shù)反演精度在一定程度上依賴于代理模型的精度，但在全局上構(gòu)建高精度的代理模型，需要通過大量的有限元正分析獲取足夠多的樣本點(diǎn)，這往往需要極大的計(jì)算消耗。若只在待反演參數(shù)點(diǎn)的附近布置較密的樣本點(diǎn)，而在其他地方稀疏取樣，使得代理模型僅在局部具有相對(duì)高的精度，便能夠通過較少的樣本點(diǎn)獲得盡可能高的反演精度。

粒子群算法PSO(Particle Swarm Optimization)是一種基于群體的隨機(jī)優(yōu)化算法[13，14]，具有較好的全局尋優(yōu)能力，但因?yàn)槠潆S機(jī)性，在局部的收斂速度較慢。基于梯度信息的高斯-牛頓法雖然容易陷入局部極小，對(duì)初值的選取要求較高，但其在局部區(qū)域往往比粒子群等啟發(fā)式算法具有更快的收斂速度以及求解精度。因此將兩者結(jié)合起來，便能夠在一定程度上同時(shí)保證全局最優(yōu)性、局部精度和收斂效率。

對(duì)于大多數(shù)參數(shù)反演問題，待反演參數(shù)總是落在參數(shù)空間的局部區(qū)域，而在全域構(gòu)建高精度的代理模型則需要消耗大量的計(jì)算機(jī)時(shí)。為此本節(jié)構(gòu)建了一種基于代理模型的迭代更新反演算法，能夠在提高反演效率的同時(shí)獲得更高的精度。

迭代更新反演的步驟如下。

(1) 在全域稀疏采樣以盡可能廣地覆蓋參數(shù)空間，通過有限元等數(shù)值方法獲得相應(yīng)的位移響應(yīng)，得到初始樣本集S0。

(2) 根據(jù)樣本集S0建立代理模型。此時(shí)的代理模型雖然精度不高但能反映總體特征。

(3) 基于PSO算法的全局尋優(yōu)能力找到一個(gè)極值點(diǎn)x0作為反演結(jié)果。由于代理模型的精度較低，反演出的結(jié)果通常和真實(shí)值也有較大的差距。

(4) 計(jì)算反演參數(shù)下的位移響應(yīng)，迭代更新。若計(jì)算出的位移與實(shí)測位移差距較大，則將此反演參數(shù)和有限元計(jì)算的位移響應(yīng)作為新的樣本對(duì)加入到初始樣本集S0中，得到新的樣本集S1。

(5) 使用增加樣本點(diǎn)后的樣本集S1重新構(gòu)建代理模型。此模型將具有在x0附近更高的精度。由于已通過PSO算法找到全局最優(yōu)的局部區(qū)域，此時(shí)可以改用局部收斂速度更快的高斯-牛頓法繼續(xù)反演，這樣能夠在提高反演速度和精度的同時(shí)避免PSO算法每次都進(jìn)行全局搜索帶來的不穩(wěn)定性。

基于POD -RBF代理模型結(jié)合PSO算法和高斯-牛頓法的迭代更新反演算法具體流程如圖1所示。

圖1 迭代更新反演算法流程

4 算例分析

4.1 模型概況

圖2 混凝土大壩有限元模型

4.2 有限元代理模型

代理模型的精度一定程度上取決于樣本數(shù)量，充足的樣本集有利于構(gòu)造高精度的代理模型，但相應(yīng)地也會(huì)消耗更多機(jī)時(shí)。因此，選取合適的樣本數(shù)量以更經(jīng)濟(jì)的方式構(gòu)建代理模型是十分重要的。

本節(jié)分別選取了不同數(shù)量的樣本集來構(gòu)建代理模型，并統(tǒng)計(jì)了代理模型對(duì)測試樣本的計(jì)算位移和對(duì)應(yīng)的有限元計(jì)算位移之間的平均相對(duì)誤差，結(jié)果如圖3所示?？梢钥闯?，針對(duì)本算例的有限元模型，當(dāng)選取樣本數(shù)量小于9時(shí)，樣本數(shù)量的增加對(duì)代理模型精度的提升作用明顯，而當(dāng)樣本數(shù)量大于9之后，樣本數(shù)量的進(jìn)一步增加并不能再顯著提升代理模型的精度。因此本算例選取9個(gè)樣本作為初始樣本集。

圖3 不同樣本個(gè)數(shù)下POD -RBF模型的平均相對(duì)誤差

在對(duì)位移響應(yīng)矩陣進(jìn)行奇異值分解時(shí)，其前幾階奇異向量包含了絕大多數(shù)特征信息，通常只需提取前幾階就足以獲得較高精度。在本算例中，對(duì)9個(gè)樣本組成的位移矩陣進(jìn)行分解，得到其前五個(gè)奇異值列入表1。

表1 位移響應(yīng)矩陣奇異值

由表1可知，位移矩陣前3階奇異值占總奇異值比重已達(dá)95%以上。截取不同個(gè)數(shù)奇異向量作為基向量進(jìn)行插值引起的平均相對(duì)誤差如 圖4 所示。代理模型的平均相對(duì)誤差的變化趨勢與奇異值的變化趨勢相反，這表明隨著奇異值占比的增加，模型精度相應(yīng)提高，當(dāng)使用前3階奇異向量構(gòu)建代理模型時(shí)，代理模型相對(duì)誤差已降至接近最低水平。

4.3 迭代更新反演結(jié)果

圖4 平均相對(duì)誤差和奇異值占比隨截取奇異值個(gè)數(shù)的變化

為驗(yàn)證本文反演方法的高效性，將文獻(xiàn)[15]的基于BP神經(jīng)網(wǎng)絡(luò)的直接反演方法加入進(jìn)行對(duì)比。本文BP神經(jīng)網(wǎng)絡(luò)的隱含層神經(jīng)元個(gè)數(shù)取為10，激活函數(shù)為Sigmoid函數(shù)。兩種方法的反演結(jié)果列入表3，由表3可知，與BP直接反演方法相比，基于POD -RBF代理模型的迭代更新反演方法在反演精度和穩(wěn)定性上表現(xiàn)都要更優(yōu)。

表2 迭代更新與無迭代更新反演結(jié)果對(duì)比

表3 迭代更新反演與BP直接反演結(jié)果對(duì)比

5 結(jié) 論

本文針對(duì)大型復(fù)雜結(jié)構(gòu)參數(shù)反演時(shí)通常耗費(fèi)大量正分析機(jī)時(shí)的問題，嘗試構(gòu)建本征正交分解(POD)結(jié)合徑向基函數(shù)(RBF)的代理模型，并結(jié)合粒子群算法(PSO)的全局尋優(yōu)能力和高斯-牛頓法(GNM)的快速局部收斂特性，開展兼具精度和效率的迭代更新反演算法探索，并應(yīng)用于典型混凝土重力壩壩段的分區(qū)彈性模量反演。

通過比較不同樣本數(shù)量下代理模型的精度，在確保代理模型準(zhǔn)確性的基礎(chǔ)上可選取合適的樣本數(shù)量，有效降低代理模型的計(jì)算機(jī)時(shí)消耗。

典型算例分析表明，基于本文方法，可先使用相對(duì)較少的樣本數(shù)量構(gòu)建相對(duì)粗糙的代理模型，然后在反演過程中不斷在全局最優(yōu)點(diǎn)附近加入新樣本更新代理模型，有針對(duì)性地提高代理模型的局部精度。通過結(jié)合粒子群算法的全局尋優(yōu)能力和高斯-牛頓法局部快速收斂能力，可用相對(duì)較少的迭代步獲得足夠高的反演精度。