帶有彈性邊界支撐梁的多宗量反問題數(shù)值求解

楊海天， 楊 博， 李哈汀

（大連理工大學(xué)工業(yè)裝備結(jié)構(gòu)分析國家重點(diǎn)實(shí)驗(yàn)室，遼寧大連 116024）

0 引 言

如何根據(jù)實(shí)際情況，合理考慮力學(xué)建模中邊界條件的影響，是一個(gè)值得探討的問題［1］.建模過程中邊界條件的不準(zhǔn)確，將會(huì)導(dǎo)致結(jié)果的偏差.以懸臂梁為例，實(shí)際情況往往可能是在約束端允許有小量位移與轉(zhuǎn)角（可被當(dāng)成一種彈性支撐）.如單純將約束假設(shè)為固支，則計(jì)算結(jié)果會(huì)出現(xiàn)較大偏差.這時(shí)將彈性支撐而不是固支作為邊界條件可能更為合理.問題是如何確定這些等價(jià)彈性支撐的彈性參數(shù).

求解帶有彈性元件邊界條件的反問題，是確定這些彈性參數(shù)的方法之一.對此國內(nèi)外已有一些相關(guān)的文獻(xiàn)報(bào)道，如文獻(xiàn)［2］以梁在靜力作用下的撓度為已知信息，利用神經(jīng)網(wǎng)絡(luò)對支撐彈簧剛度進(jìn)行了識(shí)別，但需要大量的網(wǎng)絡(luò)訓(xùn)練以提高精度，且沒有考慮噪聲影響；文獻(xiàn)［3］提出以模態(tài)參數(shù)作為目標(biāo)函數(shù)用牛頓法識(shí)別邊界支撐剛度和邊界質(zhì)量，但沒有考慮彎曲相關(guān)的彈性支撐；文獻(xiàn)［4］在支撐剛度的識(shí)別中，對牛頓法和遺傳算法進(jìn)行了比較，發(fā)現(xiàn)牛頓法對初始解的選取較敏感；也有不少學(xué)者用頻率響應(yīng)函數(shù)法識(shí)別連接參數(shù)［1］，早期的頻率響應(yīng)函數(shù)（FRF）法對連接處的噪聲較敏感，文獻(xiàn)［5］由此進(jìn)行了改進(jìn)；文獻(xiàn)［6］探討了坐標(biāo)縮減對頻率響應(yīng)函數(shù)法識(shí)別結(jié)果的影響；文獻(xiàn)［7］對動(dòng)力彈性邊界優(yōu)化設(shè)計(jì)的有關(guān)問題進(jìn)行了綜述.

以上工作主要圍繞梁彈性邊界上相關(guān)參數(shù)的反演展開，似未見到彈性邊界相關(guān)參數(shù)與梁本構(gòu)參數(shù)組合識(shí)別的文獻(xiàn)報(bào)道，而實(shí)際問題中，梁的本構(gòu)參數(shù)可能是未知的或在使用過程中發(fā)生變化.本文通過帶有彈性邊界梁的靜力有限元方程，推導(dǎo)出梁的彈性邊界／本構(gòu)參數(shù)多宗量靜力反演數(shù)值模型，采用L－M（Levenberg－Marquardt）算法［8～12］進(jìn)行求解，可同時(shí)對彈性邊界條件及梁的本構(gòu)參數(shù)進(jìn)行識(shí)別.

1 控制方程及其有限元離散

圖1所示為Bernoulli－Euler梁，梁長為L，集中力f作用于x＝L0處.

圖1 帶有彈性支撐的梁Fig.1 A beam with elastic supports

其控制方程和邊界條件為［13］

其中E為彈性模量，I為慣性矩，u為梁的橫向位移，k1和k2為橫向彈性支撐剛度，g1和g2為彎曲彈性支撐剛度.

系統(tǒng)的有限元方程可寫為［14］

其中K、Y和F分別代表總剛度陣、總節(jié)點(diǎn)位移向量和總節(jié)點(diǎn)荷載向量.

其中G為單元節(jié)點(diǎn)自由度轉(zhuǎn)換矩陣［13］，l為單元長度，N為二節(jié)點(diǎn)Hermite單元的形函數(shù)矩陣.

2 模型求解

未知參量ki、gj、Ek的識(shí)別可通過以下泛函的極小化實(shí)現(xiàn).

式中：R（X）＝y(tǒng)p（X）－ym；yp（X）、ym分別為節(jié)點(diǎn)位移向量的預(yù)測值和測量值；ki和g j為未知彈性支撐參數(shù)；Ek為任意未知單元的彈性模量.

計(jì)算流程為［8］

Step 1 選取變量初值X0，設(shè)置迭代步n＝0及ε（誤差容限）；

Step 2 計(jì)算R（X n）及R（X n），通過式（12）求解D n；

Step 3X n＋1＝X n＋D n，并令n＝n＋1.

Step 4 如果‖D n‖≤ε，則停止迭代；否則執(zhí)行Step 2.

3 噪聲處理

按照文獻(xiàn)［15］的方法注入噪聲：

其中為不含噪聲的已知信息，ξ是一個(gè)服從標(biāo)準(zhǔn)正態(tài)分布的隨機(jī)變量，Ns表示噪聲水平.對于每一個(gè)給定的Ns都由ξ生成40組y＊.

置信區(qū)間為其中表示反演結(jié)果的平均值，S表示標(biāo)準(zhǔn)差，t表示N－1個(gè)自由度的t－分布，N表示樣本容量，置信度為1－β.

4 數(shù)值算例

算例1 邊界條件的影響

圖2 彈性支撐懸臂梁Fig.2 A cantilever beam with elastic supports

考慮圖2所示帶有彈性支撐的懸臂梁，梁長L＝0.5 m，單元長度l＝0.01 m，彈性模量E＝210 GPa，慣性矩I＝10－6m4，密度ρ＝7 860 kg／m3，荷載大小f＝5 000 N.u和θ為帶有彈性邊界條件梁的位移和轉(zhuǎn)角，u0和θ0為相同荷載條件下固支邊界梁的位移和轉(zhuǎn)角.k為106～1010N／m，g為104～109N·m／rad［16］.表1給出了考慮不同邊界條件下的解.不難看出，隨梁端剛度的減小，邊界條件為彈性支撐的解與固支條件的解相比發(fā)生了明顯變化.

表1 位移轉(zhuǎn)角比較Tab.1 The contrast of displacement and angular displacement

算例2 彈性支撐與彈性模量的組合識(shí)別

考慮圖1所示兩端帶有彈性支撐的梁，將梁均勻劃分為50個(gè)單元，對彈性支撐的剛度k1、g1和梁的彈性模量E20進(jìn)行組合反演.Ei＝210 GPa（i＝1，2，…，19，21，…，50），k2＝1×105N／m，g2＝3×103N·m／rad，L0＝0.3 m，其他參數(shù)與算例1中相同.

采用3個(gè)點(diǎn)的轉(zhuǎn)角作為已知信息，測點(diǎn)位置和ξ的最大值見表2，計(jì)算結(jié)果如表3、4所示.

表2 測點(diǎn)位置和ξ的最大值Tab.2 The location of measuring points and the maximums ofξ

表3 10%噪聲下彈性支撐和彈性模量的組合識(shí)別Tab.3 Combined identification of elastic supports and modulus of elasticity at 10%of noise

表4 5%噪聲下彈性支撐和彈性模量的組合識(shí)別Tab.4 Combined identification of elastic supports and modulus of elasticity at 5%of noise

計(jì)算結(jié)果表明：當(dāng)按式（17）考慮噪聲的影響時(shí)，本文所提算法可有效地對彈性支撐的剛度和梁的本構(gòu)參數(shù)進(jìn)行識(shí)別，測點(diǎn)越接近待識(shí)別參數(shù)所在單元的位置識(shí)別效果越好.

5 結(jié) 語

本文的主要貢獻(xiàn)是建立了一個(gè)梁的彈性邊界／本構(gòu)參數(shù)的靜力組合反演數(shù)值模型，可同時(shí)對彈性邊界條件及梁的本構(gòu)參數(shù)進(jìn)行識(shí)別.數(shù)值驗(yàn)證表明，本文所提出的模型與方法是正確可行的.在本文工作的基礎(chǔ)上，作者正結(jié)合有關(guān)實(shí)驗(yàn)，開展梁的彈性邊界／本構(gòu)參數(shù)的動(dòng)力組合反演研究，以期為更合理地確定邊界條件／本構(gòu)參數(shù)提供有效的新方法.

［1］KIM T R，WU S M，EHMANN K F.Identification of joint parameters for a taper joint［J］.Journal of Engineering for Industry，1989，111（3）：282－287

［2］郭國會(huì)，易建偉.神經(jīng)網(wǎng)絡(luò)在梁的邊界條件識(shí)別中的應(yīng)用［J］.湖南大學(xué)學(xué)報(bào)（自然科學(xué)版），1998，25（4）：87－112

［3］劉玉明，易建偉.一種識(shí)別結(jié)構(gòu)系統(tǒng)邊界條件的新方法［J］.湖南大學(xué)學(xué)報(bào)（自然科學(xué)版），1999，26（2）：99－103

［4］易建偉，劉 霞.混凝土梁板類構(gòu)件邊界條件識(shí)別與研究［J］.湖南大學(xué)學(xué)報(bào)（自然科學(xué)版），2000，27（4）：81－87

［5］REN Y，BEARDS C F.Identification of joint properties of a structure using FRF data［J］.Journal of Sound and Vibration，1995，186（4）：567－587

［6］DAMJAN C，MIHA B.The influence of the coordinate reduction on the identification of the joint dynamic properties［J］.Mechanical Systems and Signal Processing，2009，23（4）：1260－1271

［7］高躍飛.結(jié)構(gòu)動(dòng)力學(xué)邊界條件優(yōu)化設(shè)計(jì)與工程實(shí)現(xiàn)方法研究［D］.西安：西北工業(yè)大學(xué)，2005

［8］王彥飛.反演問題的計(jì)算方法及其應(yīng)用［M］.北京：高等教育出版社，2007

［9］LEVENBERG K.A method for the solution of certain non－linear problems in least squares［J］.Quarterly Journal of Applied Mathematics，1944，11（2）：164－166

［10］MARQUARDT D W.An algorithm for leastsquares estimation of nonlinear inequalities［J］.SIAM Journal on Applied Mathematics，1963，11（2）：431－441

［11］MORE J J.The Levenberg－Marquardt algorithm：implementation and theory［C］／／WATSON G A，ed.Numerical Analysis Proceedings of the Biennial Conference Held at Dundee.Berlin：Springer－Verlag，1978：105－116

［12］LIU G R，HAN X.Computational Inverse Techniques in Nondestructive Evaluation［M］.London：CRC Press Inc.，2003

［13］張 雄，王天舒.計(jì)算動(dòng)力學(xué)［M］.北京：清華大學(xué)出版社，2007

［14］TABARROK B，CLEGHORN W L.Application of principle of least action to beam problems［J］.Acta Mechanica，2000，142（1）：235－243

［15］王登剛，劉迎曦，李守巨.二維穩(wěn)態(tài)導(dǎo)熱反問題的正則化解法［J］.吉林大學(xué)自然科學(xué)學(xué)報(bào)，2000（2）：55－60

［16］AGAPIOU J S.A methodology to measure joint stiffness parameters for toolholder－spindle interfaces［J］.Journal of Manufacturing Systems，2005，24（1）：13－20