人工邊界轉換方法

1、靜-動力分析中人工邊界轉換方法的研究摘要：通過將粘彈性動力人工邊界應用于同時考慮靜力效應和動力效應的工程算例，闡明了此類問題靜-動力分析人工邊界轉換時保證模型為靜力平衡體的必要性。通過將粘彈性靜-動力統(tǒng)一人工邊界應用于半無限空間體有限元模型的靜力分析中，驗證了靜力計算中的誤差將使模型動力分析的穩(wěn)態(tài)反應出現(xiàn)相近的誤差。在此基礎上，系統(tǒng)闡述了適用于同時考慮靜力效應和動力效應的工程問題的靜-動力分析人工邊界轉換方法。關鍵詞：人工邊界，靜力分析，動力分析，邊界轉換Abstract: Though the application of dynamic viscous-spring artificial

2、boundary to an engineering case with a consideration of both static and dynamic effect, and the application of the unified viscous-spring boundary for static and dynamic analysis to static analysis of a finite modal of half space, the problems of the applications of viscous-spring artificial boundar

3、y to this kind of engineering calculation was pointed out, and its corresponding solving method was proposed. On the base, a systematic switching method of these artificial boundaries was specified. Keywords: artificial boundary, static analysis, dynamic analysis, switching of boundaries1 前言人工邊界從廣義上

4、可分為靜力人工邊界和動力人工邊界。靜力人工邊界由來已久，通常有固定邊界、滾軸邊界等。動力人工邊界經(jīng)過幾十年的研究發(fā)展，已形成具有全局人工邊界和局部人工邊界的兩大類別，并應用于各自適應的工程計算中1。動力人工邊界發(fā)展到現(xiàn)在已有透射邊界、粘性邊界、粘彈性邊界等幾種類型。1994年，Deeks提出粘彈性人工邊界11。1998年，劉晶波等人發(fā)展了二維的黏彈性人工邊界3，又于2005年將其發(fā)展為三維時域黏彈性人工邊界4。2006年，劉晶波等人再將二維黏彈性邊界發(fā)展成一致粘彈性人工邊界及其對應的粘彈性邊界單元5，并于2007年推導了三維一致粘彈性人工邊界及等效粘彈性邊界單元6。目前對靜-動力分析的普遍做法

5、是采用靜力人工邊界和動力人工邊界分別對靜力問題和動力問題進行計算，將計算結果進行疊加后得到完整的結果1。但由于疊加原理僅在線彈性小變形范圍內適用，原則上不能應用于涉及非線性或大變形問題的分析。目前對涉及非線性或大變形問題的靜-動力分析，常用的人工邊界轉換方法主要有以下幾種：（1）靜力分析和動力分析都采用滾軸邊界或固定邊界；（2）靜力分析中采用滾軸邊界或固定邊界，動力分析采用粘彈性邊界、透射邊界、粘性邊界等人工邊界；（3）靜力分析和動力分析都采用靜-動力統(tǒng)一邊界，如粘彈性靜-動力統(tǒng)一人工邊界。對第（1）種方法，由于固定邊界使波動全部反射，已有許多文獻證明其具有放大振動效應的作用，目前已經(jīng)使用得不

6、多。劉晶波等人基于黏彈性動力人工邊界和半無限空間中靜力問題的基本解，建立了對動力問題和靜力問題均適用的三維黏彈性靜-動力統(tǒng)一人工邊界，從而上述第（3）種方法得以解決1。然而，在使用人工邊界對地下結構進行動力分析時，還存在一些問題。如第（2）種方法，由于在靜-動力分析的人工邊界轉換時的方法存在問題，致使產(chǎn)生錯誤的結果。在第（3）種方法中，將粘彈性靜-動力統(tǒng)一人工邊界應用于地下結構的靜力分析時，其解與準確值存在誤差。本文將就此兩問題進行論證和分析，并闡述合理的地下結構靜-動力分析人工邊界轉換方法。2 靜力和動力有限元分析原理2.1 靜力分析原理地鐵等地下工程初始應力場的確定須先計算未開挖狀態(tài)下圍巖

7、的自重應力場，進而根據(jù)施工步驟，采用釋放荷載法，計算出襯砌結構和圍巖的靜應力場13。許多地下結構的自重應力場模型可以假設為半無限空間體，根據(jù)經(jīng)典圍巖壓力理論和彈性力學理論，半無限空間體中距地表面任一深度h處的應力狀態(tài)可定義為78： (1) (2)式中，V 為豎向應力，為圍巖重度，H為橫向應力，為側壓力系數(shù)。對淺層圍巖，可假設其為各向同性介質，側壓力系數(shù)可用泊松比表示如下7： (3)對有限元計算中所取的有限區(qū)域，可以據(jù)此確定有限域邊界條件。2.2靜力分析邊界條件根據(jù)上述靜力分析方法，在靜-動力共同作用問題的計算中，在進行動力分析之前須先確定地下結構模型的靜應力場。地下結構開挖前可將大地假設為半無

8、限空間體，其在重力作用下的靜力計算，根據(jù)對稱性，模型中任一處的水平位移。因此，計算該應力場時，有限區(qū)域模型兩側可用水平約束即法向約束，底部可用全約束或僅約束豎直方向，頂面即地面應為自由邊界。2.3 動力分析原理在動力荷載作用下，有限元體系在t+Dt時刻的運動平衡方程為： (4)式中為體系的總質量矩陣；為體系的總阻尼矩陣；為體系的總剛度矩陣；為體系的節(jié)點加速度向量；為體系的節(jié)點速度向量；為體系的節(jié)點位移向量；為外荷載向量。體系的總阻尼矩陣采用瑞利阻尼9： (5)式中a、b為常數(shù)，可按兩種不同的振動頻率下測得的阻尼比加以確定。計算中常數(shù)a、b可由+w2i=2wixi和+w2j=2wjxj求得9。則

9、a和b可表示為： (6) (7)wi 、xi 分別為振型向量fi對應的自振圓頻率和阻尼比。根據(jù)振型分析結果可求得wi和wj，阻尼比在計算中取xi =xj=0.05。2.4 粘彈性人工邊界粘彈性人工邊界從用途上可分為動力人工邊界和靜-動力統(tǒng)一人工邊界；從具體實現(xiàn)方法上可分為彈簧-阻尼器邊界單元和一致粘彈性邊界單元。下面先介紹動力人工邊界。具體計算中模型邊界材料參數(shù)由其相鄰的圍巖介質材料決定。則當人工邊界采用等效的彈簧和阻尼器物理元件來模擬時，其彈簧系數(shù)和阻尼系數(shù)的計算分別如下34：法向邊界：， (8)切向邊界：， (9)式中KT、KN分別為法向與切向的彈簧剛度；R為波源至人工邊界點的距離；cs和

10、cp分別為S波和P波波速；E和G分別為介質彈性模量和剪切模量；為介質泊松比；為介質質量密度；T與N分別為切向與法向粘彈性人工邊界參數(shù)，具體取值情況見表1。表1 粘彈性動力人工邊界中參數(shù)的取值模型類型方向二維人工邊界平面內法向2.0平面內切向1.5出平面切向0.5三維人工邊界法向4.0切向2.0若人工邊界采用一致粘彈性邊界單元來模擬，其邊界單元等效剪切模量、等效彈性模量和阻尼系數(shù)可分別用如下幾式計算56： (10) (11) (12)其中等效阻尼系數(shù)取的是法向和切向的平均值。式中h為等效邊界單元厚度，為等效泊松比，n為計算模型維數(shù)。其余符號意義同上。表2 粘彈性靜-動力統(tǒng)一人工邊界中參數(shù)的取值模

11、型類型人工邊界位置方向二維人工邊界模型底面法向*/3切向2/3模型側面法向*/4切向1/8三維人工邊界模型底面法向*切向2模型側面法向*切向1/2考慮到實施的方便，實際計算中邊界材料常采用各向同性材料，這時上述等效剪切模量和等效彈性模量之間存在隱含關系式=2（1+）。考慮到普通有限元材料泊松比應限制在00.5范圍內，則該等效泊松比可按如下取值5。 (13)即計算等效泊松比時應先確定大小范圍，再確定其計算式。靜-動力統(tǒng)一人工邊界是在上述動力人工邊界的基礎上對人工邊界參數(shù)進行調整，其具體取值見表21。其中參數(shù)*具體計算公式： (10)式中為圍巖泊松比，d為位置坐標，R為荷載作用點到人工邊界點的距離

12、。對于底面人工邊界，式中D取荷載作用點至邊界單元的水平距離；對于側面人工邊界，式中D取荷載作用點至邊界單元的垂直距離。3 靜-動力分析中的人工邊界轉換的幾種方法在引言中已經(jīng)介紹了地下結構靜-動力分析中常用的三種人工邊界轉換方法，下面采用算例對第二種和第三種方法存在的問題進行分析。3.1算例一某雙線鐵路隧道，其凈空跨度10.2m，為三心圓斷面；襯砌初期支護15cm，二次襯砌30cm；隧道埋深45m，類圍巖。采用ANSYS大型有限元分析軟件對該隧道進行建模和計算，模型區(qū)域豎向取100m，橫向取100m。隧道及圍巖物理力學參數(shù)見表4。靜力分析采用的邊界條件如圖1所示，其中g為重力加速度。動力分析采用

13、如圖2如示動力模型，其最外一層即為粘彈性動力人工邊界，厚度2米。圖中P(t)為某實測的列車振動荷載，其時程曲線見圖3。 第（2）方法在靜力分析中采用滾軸邊界或固定邊界（可稱之為傳統(tǒng)靜力邊界條件），在動力分析中采用粘彈性動力人工邊界，在靜力分析和動力分析中的荷載和邊界條件見表3。 gP(t) g 圖1隧道靜力分析邊界條件 圖2 動力分析邊界條件表3 計算步驟及條件表分析步驟1靜力分析2動力分析邊界條件兩側采用滾軸邊界，底部采用固定邊界如圖1所示，一致粘彈性動力人工邊界。兩側邊界單元外邊緣法向約束，底部邊界單元外邊緣固定。荷載條件重力加速度。重力加速度g。初始應力場（靜力計算結果）。動力荷載P(t

14、)。表4 圍巖物理力學參數(shù)表材料彈性模量E泊松比粘聚力c內摩擦角密度阻尼系數(shù)圍巖3.2×109Pa0.324.5×105Pa33°2150kg/m3初期支護27.5×109Pa0.22500 kg/m3二次襯砌30×109Pa0.22500 kg/m3動力邊界9.7×107Pa00.0546圖3 列車荷載時程曲線由于隧道結構是對稱的，整個模型有限域以及荷載也是對稱的，監(jiān)測點可以取自其對稱軸及其一側即可；取隧道襯砌的拱頂、拱腰、拱腳以及仰拱底等幾個點監(jiān)測其動力響應情況。各監(jiān)測點位移時程曲線見圖4，第一主應力時程曲線見圖5。從圖4和圖5可

15、以看出，各監(jiān)測點的豎向位移和第一主應力均呈現(xiàn)首先大幅波動，然后逐漸趨于穩(wěn)定的情形。而從圖2可以看出，列車振動荷載在整個時間段為振幅相近的周期荷載，且其大小不足以使結構應力和位移發(fā)生如此大的波動。據(jù)此，可以推斷模型在零時刻就受到一個或一組與自重相當?shù)牟黄胶饬?，追溯前面的計算過程，可知這不平衡力就是靜力計算結果中的約束反力。這說明這種方法存在很大問題。圖4 算法1各監(jiān)測點豎向位移時程曲線（a）拱頂（b）拱腰（c）拱腳（d）仰拱底圖5 算法1各監(jiān)測點第一主應力時程曲線3.2 方法三的算例下面將粘彈性靜-動力統(tǒng)一人工邊界應用于方法二算例中所述的隧道模型動力分析，檢驗方法三的適用性。隧道模型如圖6所示，

16、具體的計算條件見表5，其中動力荷載P(t)為與方法二的算例相同的列車振動荷載。監(jiān)測點的設置同該算例，計算結果見圖7和圖8。表5 方法三計算步驟和條件表分析步驟1靜力分析2動力分析邊界條件施加粘彈性靜-動力統(tǒng)一人工邊界，兩側邊界單元外邊緣法向約束，底部邊界單元外邊緣固定。邊界同靜力分析邊界條件。荷載條件重力加速度g。重力加速度g。初始應力場（靜力計算結果）。動力荷載P(t)。 gP(t)圖6 雙線鐵路隧道統(tǒng)一人工邊界模型圖7 模型幾何中心點橫向應力和豎向應力圖圖8 模型底部中點橫向應力和豎向應力圖從圖中可以看出，對雙線鐵路隧道施加粘彈性靜-動力統(tǒng)一人工邊界進行動力分析時，模型節(jié)點的第一主應力和豎

17、向位移在振動初期的波動仍為整個過程的最大波動，但其波動幅值與振動荷載引起的穩(wěn)態(tài)反應幅值相近或僅有幾倍之差。這說明該方法中從靜力分析到動力分析的應力場和荷載等過渡得較好。3.3 粘彈性靜-動力統(tǒng)一人工邊界在靜力分析中的應用方法三在靜-動力分析中人工邊界轉換時靜應力場和荷載等過渡得較好，但靜力分析時沿須考慮粘彈性靜-動力統(tǒng)一人工邊界的精度，以下通過某靜力分析算例檢驗其精度。由于半無限空間體具有經(jīng)典彈性力學解，即使在巖土力學中也有理論解，前面已有論述；因此，對半無限空間體的有限元模型進行靜力分析，將結果與理論解進行對比，可以檢驗其計算方法的有效性。該模型的理論解可以通過(1)式、(2)式和(3)式算

18、得。為了比較，進行了理論解、固定邊界的數(shù)值解和粘彈性靜-動力統(tǒng)一人工邊界的數(shù)值解三種工況的計算，其計算模型分別如圖9、圖10和圖11所示。具體的計算步驟和計算條件見表6。在模型對稱軸上設置觀測點，由地面向下每10m設置一個觀測點，并將三種模型計算結果進行比較，比較結果見圖12和圖13。表6 三種工況靜力計算條件表工況理論解傳統(tǒng)靜力邊界解粘彈性靜-動力統(tǒng)一人工邊界解邊界條件如圖X所示，地面為自由邊界，地面以下為無窮遠邊界如圖X所示，兩側滾軸邊界，底部固定邊界。施加粘彈性人工邊界，兩側邊界單元外邊緣法向約束，底部邊界單元外邊緣固定。荷載條件重力加速度。重力加速度。重力加速度。 V=0H g 圖9

19、半無限空間體理論解模型（工況一） g 圖10 傳統(tǒng)靜力邊界模型（工況二） g圖11 粘彈性靜-動力統(tǒng)一人工邊界模型（工況三）圖12 豎向應力對比圖圖13 橫向應力對比圖從圖中可以看出，該有限元模型在施加粘彈性靜-動力統(tǒng)一人工邊界時，其靜力計算結果與理論解相比有較大誤差，特別是橫向應力結果誤差較大；而其施加前面所述靜力邊界條件時，豎向應力和橫向應力結果都與理論解符合得相對較好。這說明粘彈性靜-動力統(tǒng)一人工邊界應用于受重力作用的半無限空間體的靜力計算時，其誤差偏大；此結果可推廣到許多地下結構的靜力分析。4 方法二的改進根據(jù)以上分析，上述三種靜-動力分析中的人工邊界轉換方法，方法三雖然從靜力分析到動

20、力分析的應力場和荷載等過渡得較好，但其相應人工邊界靜力分析誤差偏大，本文將方法二算例中的計算程序加以改進， 將表3中動力分析的荷載條件再加上靜力分析中的約束反力。靜力計算模型仍采用圖1所示模型，動力計算采用圖14所示模型。監(jiān)測點設置與前面相同，各監(jiān)測點動力響應見圖15和圖16。 gP(t)圖14 方法二的改進動力分析模型圖15 改進后動力分析模型的監(jiān)測點豎向位移響應(a)拱頂?shù)谝恢鲬憫?b)拱腰第一主應力響應(c)拱腳第一主應力響應(d)仰拱底第一主應力響應圖16 改進后動力分析模型的應力響應從圖中可以看到，各監(jiān)測點的最大豎向位移與該方法改進前的結果相比大為減小，其后其的穩(wěn)態(tài)振動與初期相對

21、較大的振動差值比前面也小許多，第一主應力也有類似的現(xiàn)象。初期出現(xiàn)相對較大的波動是因為這里的列車振動荷載是作為突加荷載施加于模型的，上述動力模型改進后其初期動力響應波動值與穩(wěn)態(tài)反應值接近，是人工邊界轉換更為合理的表現(xiàn)。5 幾種方法的比較從以上算例中可以看到，幾種方法下各監(jiān)測點在振動初期的振幅都是最大的，這是因為這些計算中列車振動荷載都是在突加荷載的形式施加于模型的。但各種方法下振動幅大小有所不同，表7列出了各種方法下監(jiān)測點的第一主應力振幅大小。從表中可以看到，使用方法二時多數(shù)監(jiān)測點應力的振幅與使用另外兩種方法時的結果比較為最大，其主要原因是方法二在靜-動力分析的人工邊界轉換中未能較好地將靜力狀態(tài)

22、施加于模型。因此，方法二存在缺陷。表7 各種方法下的第一主應力振幅項目方法二方法三改進的方法拱頂37397.85430130拱腰102600689762拱腳59113087708930仰拱底39164.34317011375.9表8 各種方法下第一主應力的振動穩(wěn)定值與靜力解監(jiān)測點拱頂拱腰拱腳仰拱底方法二-70434.4-301746-1717520-49078.6方法三3250880-324200-20947803216940改進的方法-104182-295537-1743830-89789.9靜力解-104100-295972-1749080-82949.2從方法三的算例以及粘彈性靜-動力統(tǒng)

23、一人工邊界在靜力分析中的應用可以看到，方法三雖然從靜力分析到動力分析的應力場和荷載等過渡得較好，但其相應人工邊界應用于靜力分析時誤差偏大。表8為各種方法下各監(jiān)測點下第一主應力的振動穩(wěn)定值與靜力解的對比，此處的振動穩(wěn)定值為前面各監(jiān)測點第一主應力響應圖中呈水平線的時間段的應力值，取最后一段穩(wěn)定值。從中可看到方法三算出的應力穩(wěn)定值與靜力解差值為幾種方法中最大的，這再次證實了上述結論。從方法二改進后的算例中可以看到，該方法改進后從靜力分析到動力分析的應力場和荷載等過渡更加合理。而表8顯示改進后的方法其應力振動穩(wěn)定值與靜力解相對較為接近。因此，本文推薦改進后的方法為地下結構靜-動力分析方法。5 結論本文

24、通過將粘彈性動力人工邊界應用于同時考慮靜力效應和動力效應的工程問題，闡明了此類問題靜-動力分析人工邊界轉換時保證模型為靜力平衡體的必要性。通過將粘彈性靜-動力統(tǒng)一人工邊界應用于半無限空間體有限元模型的靜力分析中，驗證了該邊界條件及其相應計算方法在地下結構靜-動力分析中的局限性。通過對幾種方法下對同一隧道結構模型的動力分析結果的比較，本文推薦了改進的方法為最適用于地下結構靜-動力分析的方法。在此基礎上，本文提出適用于地下結構靜-動力分析的人工邊界轉換方法，其具體步驟如下：(1) 采用靜力方法和靜力約束條件對模型進行靜力分析，計算結果應包括應力場和約束反力。(2) 將靜力模型中的約束條件去掉，施加粘彈性人工邊界，在邊界單元外邊緣施加約束。(3) 在原靜力約束被去掉的地方施加約束反力，對模型施加靜力場和靜力荷載（如重力加速度或其他恒力）。(4) 對模型施加動力荷載，進行動力計算。該方法為此類工程問題在靜力分析階段和動力分析階段都使用精度較高的邊界條件提供了