大壩結(jié)構(gòu)靜動(dòng)力分析的精細(xì)化模型

林皋　龐林

摘要：闡述近期大連理工大學(xué)抗震研究所在大壩抗震分析方面所取得的研究成果，為了提高大壩結(jié)構(gòu)靜動(dòng)力分析的計(jì)算精度與效率，提出多邊形單元的計(jì)算技術(shù)與奇異應(yīng)力場(chǎng)計(jì)算方法，該方法使應(yīng)力計(jì)算成果較有限元法有很大的改進(jìn)，且計(jì)算工作量比較節(jié)省。

關(guān)鍵詞：靜動(dòng)力分析；多邊形單元；奇異應(yīng)力場(chǎng)

中圖分類號(hào)：TU45 文獻(xiàn)標(biāo)識(shí)碼：A 文章編號(hào)：1000-0666（2016）01-0001-09

0 前言

大壩結(jié)構(gòu)的安全對(duì)城市抗震防災(zāi)有重要影響。2008年5月12日發(fā)生的汶川8.0級(jí)特大地震對(duì)位于成都市西北60km岷江上游的紫坪鋪大壩造成了一定的地震損傷，威脅到成都市的安全。吉林市的安全與豐滿大壩的安全息息相關(guān)，密云水庫、官廳水庫的安全也對(duì)北京市的安全產(chǎn)生一定影響。研究城市防災(zāi)需要關(guān)注周邊大壩的安全設(shè)防。

新中國成立以來，特別是改革開放以來，我國水利水電建設(shè)迅速發(fā)展。進(jìn)入21世紀(jì)，世界大壩建設(shè)中心已經(jīng)轉(zhuǎn)向中國。一大批世界級(jí)的高壩已經(jīng)或正在中國大地上興建，壩高屢破世界紀(jì)錄，其成就舉世矚目。我國地處世界兩大活躍地震帶——環(huán)太平洋地震帶和歐亞地震帶的交匯地段，地震活動(dòng)頻繁。為適應(yīng)地震危險(xiǎn)區(qū)大壩建設(shè)的需要，我國的大壩抗震技術(shù)也蓬勃發(fā)展，并進(jìn)入世界先進(jìn)行列。

自20世紀(jì)50年代以來，大連理工大學(xué)工程抗震研究所一直根據(jù)我國大壩抗震的需要，致力于研究如何提高大壩抗震分析的水平。大壩抗震的特點(diǎn)在于大壩所處環(huán)境復(fù)雜（圖1），參照世界壩工專家Chopra（1992，2008）的意見，大壩抗震的關(guān)鍵技術(shù)包括壩與庫水的動(dòng)力相互作用分析，壩與無限地基的動(dòng)力相互作用分析，筑壩材料的動(dòng)態(tài)特性及大壩地震損傷發(fā)展與風(fēng)險(xiǎn)分析。在2013年中國工程院舉辦的“中國工程科技論壇——水安全與水利水電可持續(xù)發(fā)展”會(huì)議上，林皋（2014）提出了大壩抗震模型的改進(jìn)意見。關(guān)于壩與庫水的動(dòng)力相互作用分析需要考慮水的可壓縮性以及動(dòng)水壓力波在水庫邊界的吸收作用，Chopra（2008）采用有限元法（FEM）進(jìn)行求解，動(dòng)水壓力按壩的振動(dòng)模態(tài)展開，計(jì)算復(fù)雜、難以實(shí)際應(yīng)用。Lin等（2007a，2012a）采用比例邊界有限元法（SBFEM）進(jìn)行求解，按庫水的模態(tài)進(jìn)行展開，計(jì)算簡(jiǎn)便、精度高、便于實(shí)際應(yīng)用。關(guān)于壩與無限地基的動(dòng)力相互作用，由于計(jì)算模型比較復(fù)雜，文獻(xiàn)中一般都將無限地基簡(jiǎn)化為均質(zhì)地基進(jìn)行考慮，難以反映實(shí)際情況。我們將無限地基分為層狀分布、分塊不均質(zhì)和近場(chǎng)雜亂不規(guī)則等幾種類型，分別提出了不同的計(jì)算模型（Linet al，2015，2007b；Yin et al，2013），計(jì)算結(jié)果都具有較好的精度與效率。為了提高壩-庫水-地基系統(tǒng)時(shí)域分析的計(jì)算精度與效率，Lin等（2012b）提出以等幾何分析方法（isogeometric a-nalysis，簡(jiǎn)稱IGA）代替FEM來進(jìn)行壩的應(yīng)力分析，并將IGA與SBFEM相結(jié)合來進(jìn)行壩與庫水和壩與無限地基的動(dòng)力相互作用分析。關(guān)于壩材料的動(dòng)力特性方面，混凝土是速率敏感性材料。自從第一顆原子彈1945年在日本廣島爆炸以后，為了防衛(wèi)核爆炸造成的災(zāi)難，許多國家都加強(qiáng)了混凝土材料抗擊核沖擊的研究。但這種研究偏重單調(diào)瞬時(shí)沖擊，而且抗壓性能居多，抗拉性能偏少。地震荷載的特點(diǎn)是不規(guī)則的循環(huán)往復(fù)，而且加載速率也偏低。為此，我們開展了2000多試件混凝土在地震作用下的動(dòng)態(tài)特性的研究，取得了一批有價(jià)值的成果（Yan，Lin，2006；Lin et al，2007c）。

筆者認(rèn)為除了計(jì)算模型以外，數(shù)值計(jì)算方法也對(duì)大壩抗震分析的計(jì)算精度與效率產(chǎn)生影響。大壩地震響應(yīng)分析中，有限元方法的應(yīng)用比較普遍。有限元法可以適應(yīng)各種不規(guī)則的幾何邊界條件和分布不規(guī)律的單元材料特性，可以應(yīng)用于各種類型物理問題的求解，計(jì)算能力強(qiáng)，不受問題求解自由度的限制。但有限元法也有其不足之處，其主要采用低階插值函數(shù)，連續(xù)性差，應(yīng)力計(jì)算精度低。但在工程設(shè)計(jì)中應(yīng)力是對(duì)大壩安全性評(píng)價(jià)的重要依據(jù)，所以，有必要尋求更為有效的其它數(shù)值計(jì)算方法。根據(jù)我們研究，提出以下計(jì)算模型和方法供參考，歡迎批評(píng)指正。

1 多邊形單元技術(shù)

多邊形單元（PE-Polygon Element）不是一個(gè)新概念，但經(jīng)過澳大利亞新南威爾士大學(xué)（UN-SW）宋崇民教授團(tuán)隊(duì)的發(fā)展，建立在SBFEM基礎(chǔ)上的多邊形單元已經(jīng)成為與FEM一樣具有廣泛適應(yīng)性的網(wǎng)格離散技術(shù)，對(duì)數(shù)值計(jì)算具有強(qiáng)大的功能。PE可稱為廣義的FE。

大連理工大學(xué)發(fā)展了PE技術(shù)在結(jié)構(gòu)靜動(dòng)力分析中的應(yīng)用，主要工作為：（1）單元細(xì)化與連續(xù)性升階可以只在PE的邊界上進(jìn)行，不改變網(wǎng)格的拓樸結(jié)構(gòu)，使應(yīng)力計(jì)算精度顯著提高，但計(jì)算工作量則很節(jié)約；（2）將PE與曲線邊界相結(jié)合，既可構(gòu)造純直線邊的PE，也可構(gòu)造曲線邊的PE，或直線邊與曲線邊相結(jié)合的PE；（3）將PE應(yīng)用于結(jié)構(gòu)的奇異應(yīng)力場(chǎng)分析，可以達(dá)到很高的計(jì)算精度。

多邊形單元按SBFEM進(jìn)行建模（圖2）。單元域內(nèi)選擇一點(diǎn)0作為比例相似中心，其要求是從0點(diǎn)可以看見邊界上各點(diǎn)，或0點(diǎn)可與邊界上任意點(diǎn)用直線連接，沒有障礙。單元各邊上可以任意劃分子單元和結(jié)點(diǎn)，數(shù)目不限。如圖2所示，PE為六邊形，每邊可以只劃分一個(gè)子單元，每個(gè)子單元含3個(gè)結(jié)點(diǎn)，則具有2階精度。AB邊也可以劃分2個(gè)子單元，每個(gè)子單元各3個(gè)結(jié)點(diǎn)共5個(gè)結(jié)點(diǎn)，這就相當(dāng)于網(wǎng)格加密。各子單元用SBFEM進(jìn)行建模，由于SBFEM為半解析方法，從0點(diǎn)至邊界AM上任意點(diǎn)相連的射線上，變量具有解析解，這就使PE計(jì)算域具有很高的連續(xù)性，使應(yīng)力計(jì)算達(dá)到很高的計(jì)算精度。

含有曲線邊的多邊形單元如圖3所示。這時(shí)，含曲線邊AF的子單元1^e按等幾何分析（IGA）建模，采用非均勻有理B樣條（NURBS）作為插值基函數(shù)，即使在很粗網(wǎng)格的條件下也可以準(zhǔn)確地描述幾何外形。IGA充分利用了NURBS的優(yōu)點(diǎn)，可以在保持幾何外形不變的條件下方便地進(jìn)行網(wǎng)格細(xì)分，提升插值的階次，提高計(jì)算的精度，所以很適于曲線形結(jié)構(gòu)的建模。NURBS基函數(shù)由控制點(diǎn)進(jìn)行描述，不具有插值性。圖3中BE、CD、DE各邊相應(yīng)的子單元由SBFEM進(jìn)行建模，AB、EF邊相對(duì)應(yīng)的兩個(gè)子單元2^e和6^e按SBFEM建模，但是由于和1#相連需作一定調(diào)整。亦即A點(diǎn)的位移需用控制點(diǎn)1、2、3的位移加以表示，F(xiàn)點(diǎn)的位移需用控制點(diǎn)2、3、4的位移加以表示，這樣就保持了PE單元位移場(chǎng)的連續(xù)性。

PE由域內(nèi)三角形單元形心點(diǎn)連線加周邊三角形單元中點(diǎn)連線構(gòu)成（圖4），所以建模非常方便。單元便于細(xì)分和拼合，如圖5所示。

多邊形單元具有以下特點(diǎn)：（1）適于各種幾何形體的建模，具有高度靈活性；（2）單元具備單位分解性與線性完備性，通過分片檢驗(yàn)；

（3）可以將網(wǎng)格加密進(jìn)行單元細(xì)化，也可在不改變網(wǎng)格拓?fù)浣Y(jié)構(gòu)的條件下，在單元各邊上增加結(jié)點(diǎn)實(shí)現(xiàn)單元細(xì)化和連續(xù)性升階；（4）可以方便地進(jìn)行粗細(xì)網(wǎng)格過度，只需在不同網(wǎng)格交界面上增加結(jié)點(diǎn)即可（圖6）；（5）可以和FE和NuRBS兼容。

2 數(shù)值算例

為了檢驗(yàn)PE的效果，通過以下算例加以說明。

2.1 厚壁圓筒承受均勻內(nèi)壓作用

本例有解析解。利用對(duì)稱性，只取厚壁圓筒的1/4進(jìn)行離散。假設(shè)內(nèi)壓力p=1，徑向與環(huán)向應(yīng)力只和內(nèi)徑和外徑的比值有關(guān)，假設(shè)r/R=1/3。將PE與FE（二階單元）的計(jì)算結(jié)果進(jìn)行對(duì)比，兩種方法的網(wǎng)格劃分如圖7所示。兩種方法計(jì)算精度與解析解相比的誤差隨網(wǎng)格結(jié)點(diǎn)數(shù)的變化的對(duì)比如圖8所示，圖中PE代表直線邊界網(wǎng)格，CPE代表內(nèi)外曲線邊界網(wǎng)格，由圖可見PE的計(jì)算優(yōu)越性是明顯的。

2.2 曲線形溢流壩受靜水壓力作用

某工程溢流壩，高H=66m，受靜水壓力作用。壩材料參數(shù)假設(shè)為：彈性模量E=36.5GPa；泊松比γ=0.167。采用PE與FE兩種方法對(duì)該工程溢流壩進(jìn)行計(jì)算對(duì)比。壩體兩種方法網(wǎng)格劃分如圖9所示，廊道部位網(wǎng)格劃分如圖10所示。圖10圖題中p表示單元離散所用的插值階數(shù)。廊道上部曲線邊界的應(yīng)力分量對(duì)比如圖11所示，廊道周邊第一主應(yīng)力分布如圖12所示。由圖可見，有限元網(wǎng)格即使在廊道角點(diǎn)O附近細(xì)分，應(yīng)力奇異性仍無法反映。

3 異應(yīng)力場(chǎng)分析

有限元方法無法處理應(yīng)力奇異性問題。SBFEM由于徑向具有解析解，所以特別有利于處理應(yīng)力奇異性問題。關(guān)鍵是須將應(yīng)力奇異點(diǎn)設(shè)置于相似中心處。但是SBFEM中，通過相似中心的邊界不能進(jìn)行離散，只能將此邊界上作用的外力或變形轉(zhuǎn)化為相應(yīng)的邊界條件，在計(jì)算方程中加以體現(xiàn)。我們通過算例來加以說明。

研究壩踵點(diǎn)在水壓力作用下的奇異應(yīng)力分布。圖13表示PE和FE兩種網(wǎng)格劃分。

圖14表示PE與FE在壩踵點(diǎn)附近的網(wǎng)格形狀和節(jié)點(diǎn)分布。圖15表示各種方法的計(jì)算結(jié)果。圖中SBFEM相似中心位于域內(nèi)，但采用不同的結(jié)點(diǎn)數(shù)進(jìn)行離散。PE為相似中心位于壩踵A點(diǎn)，但采用不同的邊界條件。由圖可見，F(xiàn)EM即使進(jìn)行網(wǎng)格細(xì)分都無法表示奇異應(yīng)力場(chǎng)的分布，采用SBFEM時(shí)，如不將相似中心置于應(yīng)力奇異點(diǎn)處也難以描述壩踵點(diǎn)附近的奇異應(yīng)力分布。采用PE方法同時(shí)將相似中心置于壩踵點(diǎn)有可能獲得接近奇異應(yīng)力場(chǎng)分布的結(jié)果。這時(shí)需通過其它方法求得AB和AC兩邊界面上的位移，以建立相應(yīng)的AB和AC邊的Dirichlet邊界條件。但處理方式不同所取得的效果也有差別。筆者采用了兩種方式處理。PE1中AB和AC兩邊界均按其它方法計(jì)算出的位移建立Dirichlet邊界條件進(jìn)行計(jì)算，而PE2中邊界AC采用已知水壓力建立Neumann邊界條件，AB邊界仍采用其它方法計(jì)算的位移建立Dirichlet邊界條件?？梢钥闯鯬E2的效果比較好，這是因?yàn)锳C邊的邊界條件比較準(zhǔn)確，同時(shí)AC邊上的水壓力也對(duì)奇異應(yīng)力場(chǎng)影響比較大的緣故。因此，為了取得良好的奇異應(yīng)力場(chǎng)的效果，建立邊界面的準(zhǔn)確邊界條件仍然很重要。

動(dòng)力問題奇異應(yīng)力場(chǎng)的分析也可采用相同的方法進(jìn)行處理。圖16～20為地震波作用下壩踵點(diǎn)奇異應(yīng)力場(chǎng)的分析結(jié)果。采用的材料參數(shù)列入表1，為了表示PE的優(yōu)越性，將PE結(jié)果與比較準(zhǔn)確的IGA的計(jì)算結(jié)果進(jìn)行比較。IGA的各級(jí)網(wǎng)格均可自動(dòng)加密，網(wǎng)格采用2階單元。圖18表示奇異應(yīng)力場(chǎng)分析結(jié)果的比較。由圖可見，在距壩踵點(diǎn)水平距離x=1.0m以外，各種網(wǎng)格的計(jì)算結(jié)果都收斂而且相符性比較好。但距離壩踵點(diǎn)x=0.5m范圍以內(nèi)，IGA1-IGA3偏差較大，IGA4結(jié)果與本文方法吻合良好，達(dá)到收斂要求。圖19為壩踵點(diǎn)的位移與加速度響應(yīng)時(shí)程比較。由于位移不存奇性，與應(yīng)力結(jié)果不同，各網(wǎng)格計(jì)算結(jié)果相差不大。

4 結(jié)語

按比例邊界有限元法（SBFEM）原理建立的多邊形單元（PE）是一種提高應(yīng)力計(jì)算精度的高效數(shù)值方法。多邊形單元可以進(jìn)一步發(fā)展為多面體單元，即將二維問題發(fā)展到三維問題，具有廣闊的發(fā)展前景。