蔣新新，鐘 紅，李云途，牛景太，鄧志平，黃紅元

(1.南昌工程學(xué)院，江西 南昌 330099；2.流域水循環(huán)模擬與調(diào)控國家重點(diǎn)實(shí)驗(yàn)室，中國水利水電科學(xué)研究院，北京 100038)

1 研究背景

我國是一個(gè)地震多發(fā)國家，西部地區(qū)建設(shè)的多座高壩工程面臨嚴(yán)峻的抗震問題，對此工程界和學(xué)術(shù)界予以高度關(guān)注[1-2]。按照現(xiàn)行抗震規(guī)范的規(guī)定，對抗震設(shè)防類別為甲類的重力壩，除按設(shè)計(jì)地震進(jìn)行抗震設(shè)計(jì)外，通常還應(yīng)論證其在最大可信地震下的災(zāi)變安全裕度，允許出現(xiàn)一定程度的損傷開裂，但須避免發(fā)生潰壩導(dǎo)致庫水失控下泄[3-4]。在強(qiáng)震作用下重力壩產(chǎn)生嚴(yán)重震損的實(shí)例不多，其中以印度Koyna重力壩最為典型。該壩在1967年Koyna地震中非溢流壩段出現(xiàn)貫穿性裂縫，促使學(xué)術(shù)界對重力壩地震開裂破壞機(jī)理開展深入研究[5]。Chopra等[6]對Koyna大壩線彈性地震響應(yīng)的分析表明壩體上部出現(xiàn)了較大拉應(yīng)力，根據(jù)拉應(yīng)力分布預(yù)測的壩體開裂區(qū)域與觀測結(jié)果基本吻合；Batta等[7]應(yīng)用線彈性斷裂力學(xué)和邊界元分別模擬了單條和多條預(yù)設(shè)裂縫情況下的Koyna大壩裂縫擴(kuò)展過程，發(fā)現(xiàn)在壩體上下游面預(yù)設(shè)的裂縫均可能擴(kuò)展形成貫通裂縫；Zhang等[8]采用彌散裂縫模型和擴(kuò)展有限元法分析了裂縫擴(kuò)展對大壩地震響應(yīng)的影響，結(jié)果表明下游折坡處最先開裂并向上游壩面擴(kuò)展直至形成貫穿裂縫。地震作用下大壩裂縫擴(kuò)展會(huì)影響壩體的整體性，降低其整體剛度，在地震往復(fù)作用下出現(xiàn)裂縫的快速發(fā)展導(dǎo)致大壩失效破壞[9-10]。因此，為準(zhǔn)確評(píng)估大壩在強(qiáng)地震作用下的響應(yīng)，發(fā)展混凝土壩動(dòng)態(tài)開裂擴(kuò)展數(shù)值模型很有必要，其中，描述混凝土開裂引起的非連續(xù)位移場和預(yù)測裂縫擴(kuò)展路徑是關(guān)鍵問題。

混凝土壩的斷裂力學(xué)模型主要有分離式裂縫模型[11]和彌散裂縫模型[12]。分離式裂縫模型通過在裂縫擴(kuò)展過程中引入新的節(jié)點(diǎn)描述裂縫軌跡，能夠直觀體現(xiàn)裂縫的擴(kuò)展過程，且可以方便地考慮縫內(nèi)水壓力的影響。這種裂縫擴(kuò)展模擬方式在有限差分法和有限元方法中應(yīng)用較多，但是由于必須通過網(wǎng)格重剖分反映裂縫的擴(kuò)展過程，導(dǎo)致前處理和計(jì)算工作量較大。為了提升分析效率，局部網(wǎng)格重剖分是有效的解決辦法[13]。彌散裂縫模型通過降低開裂單元的剛度間接體現(xiàn)裂縫的影響，其優(yōu)勢是裂縫擴(kuò)展過程中網(wǎng)格不需重剖分，代表性方法為擴(kuò)展有限元法。目前擴(kuò)展有限元法已經(jīng)集成于ABAQUS軟件平臺(tái)，借助于ABAQUS強(qiáng)大的前后處理實(shí)現(xiàn)了多種類型的裂縫擴(kuò)展模擬。但不能忽視的是，集成過程中簡化了擴(kuò)展有限元法裂縫尖端增強(qiáng)功能，且裂縫尖端需要位于單元邊界上，可能導(dǎo)致計(jì)算不準(zhǔn)確和收斂性問題[14]，從而影響裂縫擴(kuò)展模擬的精度。

比例邊界有限元法(SBFEM)[15]作為一種半解析的數(shù)值分析方法，兼具有限元法和邊界元法的優(yōu)勢，具有降低一維和徑向解析求解的特點(diǎn)[16-17]。近年來SBFEM在斷裂分析中展現(xiàn)了明顯的優(yōu)勢[18]，大為簡化裂縫尖端奇異應(yīng)力場的描述，且廣義應(yīng)力強(qiáng)度因子直接從單元應(yīng)力模態(tài)中提取、適用于均質(zhì)材料和界面材料等多種情況[19]。此外，多邊形單元具有靈活地分割重組網(wǎng)格的優(yōu)勢，相應(yīng)發(fā)展的局部網(wǎng)格重剖分技術(shù)能夠高效地描述裂縫的擴(kuò)展過程。施明光等[20]、鐘紅等[21]基于SBFEM和多邊形網(wǎng)格局部重剖分技術(shù)分析了重力壩在靜水超載作用下的斷裂過程；Li 等[22]建立了重力壩水力劈裂的SBFEM模型，著重研究了裂縫張合過程中縫內(nèi)水壓力的變化及對大壩斷裂的影響。本文在此基礎(chǔ)上將多邊形比例邊界有限元法應(yīng)用于重力壩地震斷裂分析，推導(dǎo)了相應(yīng)的動(dòng)力控制方程；提出了廣義動(dòng)態(tài)應(yīng)力強(qiáng)度因子時(shí)域計(jì)算方法，考慮裂縫面動(dòng)態(tài)接觸條件，建立了一種全自動(dòng)的重力壩分離式裂縫動(dòng)態(tài)擴(kuò)展模型；模擬了Koyna重力壩的地震斷裂過程，并對模型網(wǎng)格密度和裂縫擴(kuò)展步長等參數(shù)進(jìn)行了敏感性分析。

2 多邊形比例邊界有限元原理與動(dòng)力方程

2.1 多邊形比例邊界有限單元基本原理比例邊界有限元法是一種在域邊界離散而域內(nèi)徑向解析的數(shù)值分析方法。采用該方法的基本要求是滿足相似條件，即站在相似中心位置，域邊界的所有點(diǎn)均直接可見。若采用多邊形對求解域進(jìn)行離散，每個(gè)多邊形采用比例邊界有限元法求解[23]，則為多邊形比例邊界有限元法。如圖1(a)所示，建立含預(yù)設(shè)裂縫的重力壩多邊形比例邊界有限元模型，多邊形單元可以方便地適應(yīng)模型復(fù)雜邊界，預(yù)設(shè)裂縫由多邊形裂縫比例邊界有限單元描述。如圖1(b)所示，裂縫單元相似中心“SC”位于裂縫尖端，整個(gè)單元由邊界向相似中心等比例縮放形成，其中裂縫面由裂縫開口節(jié)點(diǎn)向相似中心縮放來描述，無須特殊處理。多邊形比例邊界單元徑向坐標(biāo)記為ξ，在相似中心ξ=0，在單元邊界處ξ=1，取值范圍為[0，1]；邊界采用兩節(jié)點(diǎn)線性單元離散，線性單元節(jié)點(diǎn)逆時(shí)針排列，環(huán)向坐標(biāo)記為η，取值范圍為[-1，1]。

圖1 重力壩多邊形比例邊界有限元模型示意

在不計(jì)體力和裂縫面荷載情況下，由虛功原理可導(dǎo)出多邊形比例邊界有限元控制方程[15，19]為：

(1)

式中：{u(ξ)}為節(jié)點(diǎn)位移函數(shù)；{q(ξ)}為節(jié)點(diǎn)力函數(shù)；[Z]為Hamilton矩陣[19]。

采用特征值分解或者Schur分解[18]均可求解控制方程(1)。以特征值分解為例，多邊形內(nèi)任一點(diǎn)的位移{u(ξ，η)}可表示為節(jié)點(diǎn)位移函數(shù){u(ξ)}沿環(huán)向邊界的插值，具體表達(dá)式為：

(2)

由位移導(dǎo)出對應(yīng)多邊形內(nèi)任一點(diǎn)的應(yīng)變和應(yīng)力，其表達(dá)式為：

(3)

(4)

2.2 剛度矩陣和質(zhì)量矩陣對于任意的多邊形比例邊界有限單元，單元?jiǎng)偠染仃嘯Ke]表達(dá)式為：

(5)

多邊形比例邊界有限單元質(zhì)量矩陣由邊界上線性單元根據(jù)單元節(jié)點(diǎn)連接組裝得到，其表達(dá)式為：

(6)

式中[m]為系數(shù)矩陣的縮略表達(dá)式，其僅取決于多邊形單元的密度和單元形狀[15]。

2.3 多邊形比例邊界有限元靜動(dòng)力方程對于由多邊形單元離散的任意結(jié)構(gòu)體系，結(jié)構(gòu)整體剛度陣[K]和整體質(zhì)量陣[M]類似于有限元法通過單元組裝形成，因此靜力平衡方程和動(dòng)力方程分別為：

[K]syggg00={F}

(7)

(8)

實(shí)際工程動(dòng)力分析需要考慮阻尼的影響，采用瑞利阻尼，阻尼矩陣通過結(jié)構(gòu)整體剛度陣和質(zhì)量陣的線性組合得到，其表達(dá)式為：

[C]=α0[M]+α1[K]

(9)

式中α0、α1為組合系數(shù)，均與結(jié)構(gòu)的頻率和材料阻尼相關(guān)。

采用Newmark方法求解式(8)可獲得任意時(shí)刻結(jié)構(gòu)的位移場、速度場和加速度場，應(yīng)變場和應(yīng)力場可分別由式(3)(4)通過位移場導(dǎo)出。

從總的流程來說，對于采用多邊形比例邊界有限單元離散的結(jié)構(gòu)動(dòng)態(tài)斷裂模擬，類似于有限元法，每一時(shí)間步先進(jìn)行結(jié)構(gòu)動(dòng)力分析，在此基礎(chǔ)上求解結(jié)構(gòu)斷裂參數(shù)并進(jìn)行動(dòng)態(tài)斷裂分析。

3 廣義動(dòng)態(tài)應(yīng)力強(qiáng)度因子及裂縫擴(kuò)展準(zhǔn)則

3.1 廣義裂縫動(dòng)態(tài)應(yīng)力強(qiáng)度因子對于無限平面內(nèi)以速度為v擴(kuò)展傳播的穩(wěn)態(tài)裂縫，其動(dòng)態(tài)應(yīng)力強(qiáng)度因子[24]計(jì)算公式為：

(10)

式中：應(yīng)力強(qiáng)度因子K上標(biāo)“eq”表示依據(jù)瞬時(shí)應(yīng)力場計(jì)算的瞬態(tài)應(yīng)力強(qiáng)度因子，上標(biāo)“dyna”代表動(dòng)態(tài)應(yīng)力強(qiáng)度因子；kⅠ(v)、kⅡ(v)與裂縫傳播速度和材料特性相關(guān)，其表達(dá)式為：

(11)

式中cd、cs和cR分別為壓縮波、剪切波和瑞利波的波速。

3.2 多邊形裂縫單元廣義動(dòng)態(tài)應(yīng)力強(qiáng)度因子求解在SBFEM中，Song等[19]提出了廣義應(yīng)力強(qiáng)度因子以描述裂縫尖端應(yīng)力奇異性，對均質(zhì)材料斷裂、界面斷裂和V型缺口等斷裂類型均適用。本文將其推廣至結(jié)構(gòu)的動(dòng)態(tài)斷裂分析，建立多邊形裂縫單元的廣義動(dòng)態(tài)應(yīng)力強(qiáng)度因子求解方法。

(12)

圖2 由SBFEM模擬的多邊形裂縫單元示意

引入裂縫特征長度L，根據(jù)關(guān)系式ξ=r/L推導(dǎo)瞬態(tài)廣義應(yīng)力強(qiáng)度因子為：

(13)

(14)

將裂縫擴(kuò)展角θc代入下式，計(jì)算t時(shí)刻對應(yīng)裂縫擴(kuò)展角的等效應(yīng)力強(qiáng)度因子K(t)：

(15)

當(dāng)?shù)刃?yīng)力強(qiáng)度因子達(dá)到材料的斷裂韌度，認(rèn)為裂縫處于臨界狀態(tài)，即將擴(kuò)展。在混凝土大壩斷裂分析中裂縫擴(kuò)展長度取為固定值Δa，每次擴(kuò)展調(diào)用局部網(wǎng)格重剖分方法以形成新的裂縫邊界，當(dāng)裂尖觸及模型邊界或者動(dòng)力分析完成均會(huì)觸發(fā)斷裂分析結(jié)束。

4 裂縫動(dòng)態(tài)擴(kuò)展模擬過程

4.1 多邊形單元生成將結(jié)構(gòu)自由剖分為三角形背景網(wǎng)格，進(jìn)而轉(zhuǎn)化為多邊形單元。如圖3(a)所示，對重力壩模型自由剖分得到三角形背景網(wǎng)格。如圖3(b)所示，在模型內(nèi)部，以三角形網(wǎng)格頂點(diǎn)為相似中心，將周圍所有三角形的形心逆時(shí)針依次相連形成多邊形；對于模型邊界上的三角形頂點(diǎn)，將周圍的三角形形心及邊界線中點(diǎn)相連構(gòu)造多邊形，新生成多邊形的形心作為相似中心；在裂尖點(diǎn)生成的多邊形裂縫單元，取裂尖點(diǎn)坐標(biāo)作為裂縫單元相似中心。最終形成含預(yù)設(shè)裂縫的重力壩多邊形比例邊界有限元模型如圖3(c)所示。

圖3 三角形網(wǎng)格轉(zhuǎn)化為多邊形網(wǎng)格示意

4.2 裂縫擴(kuò)展局部網(wǎng)格重剖分裂縫擴(kuò)展局部網(wǎng)格重剖分方法是對裂尖點(diǎn)附近局部網(wǎng)格區(qū)域進(jìn)行重剖分，每一步擴(kuò)展過程所需前處理工作量很小，只需更新重剖分區(qū)域內(nèi)的節(jié)點(diǎn)及其位移、速度和加速度，從而有效地提高斷裂分析效率。具體的網(wǎng)格重剖分過程如下：

(1)圖4(a)為裂尖點(diǎn)附近三角形背景網(wǎng)格(黑線)和生成的多邊形網(wǎng)格(藍(lán)線)，C1為擴(kuò)展前的裂尖點(diǎn)，由式(14)(15)根據(jù)裂縫應(yīng)力強(qiáng)度因子計(jì)算裂縫擴(kuò)展方向，再由所設(shè)定的裂紋擴(kuò)展步長Δa得到新的裂尖點(diǎn)C2，依據(jù)裂縫擴(kuò)展軌跡C1C2(紅線)確定所切割的三角形背景網(wǎng)格(紫線包圍區(qū)域)；

圖4 多邊形網(wǎng)格局部重剖分過程

(2)找出與被切割三角形相連的周邊三角形網(wǎng)格，連接這些三角形的外邊界，所圍成的區(qū)域即為需要重新剖分的區(qū)域(圖4(b)中的紅線包圍區(qū)域)；

(3)對圖4(b)中的紅線包圍區(qū)域重新剖分三角形網(wǎng)格，其中新生成的裂縫面成為新的三角形網(wǎng)格邊界的一部分，如圖4(c)所示。將重剖分的三角形網(wǎng)格和圖4(b)中的紅線區(qū)域外原有三角形網(wǎng)格組合形成新的三角形背景網(wǎng)格，之后再按照4.1節(jié)給出的三角形轉(zhuǎn)化為多邊形網(wǎng)格過程生成新的多邊形網(wǎng)格。

4.3 裂縫動(dòng)態(tài)接觸模擬算法大壩在地震往復(fù)作用下，若不對裂縫面施加任何約束條件，可能出現(xiàn)裂縫面相互嵌入的現(xiàn)象，即裂縫面法向相對位移出現(xiàn)負(fù)值。如圖5所示，通過在裂縫面節(jié)點(diǎn)對之間施加接觸彈簧以阻止裂縫面的嵌入，其基本思想類似于罰函數(shù)法，即將接觸約束條件當(dāng)作懲罰項(xiàng)加到系統(tǒng)的總泛函中，再對總泛函求極小。動(dòng)態(tài)斷裂分析時(shí)，一旦判斷裂縫面相互嵌入，在動(dòng)力平衡方程中加入接觸彈簧剛度矩陣，不斷調(diào)整其剛度迭代求解方程，直至滿足裂縫面接觸條件。接觸彈簧剛度矩陣表達(dá)式為：

(16)

圖5 裂縫面節(jié)點(diǎn)對施加接觸彈簧示意

5 數(shù)值算例

Koyna重力壩地震破壞過程是大壩地震分析中的經(jīng)典研究對象，本文模擬了該壩在地震作用下的斷裂過程。已有研究表明裂縫最容易從下游折坡處起裂擴(kuò)展形成上下游貫穿裂縫[7-8]，故在基于斷裂力學(xué)的模擬中通常在下游折坡處預(yù)設(shè)裂縫。為與文獻(xiàn)[8]結(jié)果進(jìn)行對比，大壩的幾何尺寸、材料參數(shù)、約束條件和荷載等均與之相同。大壩斷面如圖6所示，在壩體下游面高程66.5 m處(折坡點(diǎn))預(yù)設(shè)水平裂縫，長度為1 m。壩體混凝土材料參數(shù)為：E=31 GPa，泊松比ν=0.2，密度ρ=2643 kg/m3。動(dòng)力分析中壩體混凝土動(dòng)態(tài)彈性模量取為35.7 GPa[27]，阻尼比取為5%，動(dòng)力分析步長取為0.02 s，混凝土動(dòng)態(tài)斷裂韌度KⅠd=1.8 MN/m3/2。計(jì)算中采用平面應(yīng)變假定。

圖6 預(yù)設(shè)裂縫的Koyna重力壩典型剖面圖(單位：m)

上游水位取為96.5 m。水平向、豎向加速度幅值分別為0.49g和0.34g，圖7為地震動(dòng)加速度時(shí)程。采用Westergaard附加質(zhì)量模型考慮動(dòng)水壓力，不考慮壩體-地基動(dòng)力相互作用的影響。為阻止地震作用過程中裂縫面的嵌入，在裂縫面施加接觸彈簧。

圖7 Koyna地震波加速度時(shí)程

5.1 靜力計(jì)算結(jié)果靜力分析為后續(xù)動(dòng)力分析提供初始狀態(tài)。靜荷載包括壩體自重和上游靜水壓力，分別采用圖8所示的粗網(wǎng)格和細(xì)網(wǎng)格方案計(jì)算大壩響應(yīng)，兩種網(wǎng)格方案(圖中G為壩體自重，左邊箭頭表示水壓力)分別包含190個(gè)和475個(gè)單元，以及439個(gè)和1043個(gè)節(jié)點(diǎn)。圖9為大壩變形圖，其中圖9(a)(b)均對應(yīng)設(shè)置縫面接觸彈簧的結(jié)果，圖9(c)未設(shè)置縫面接觸彈簧。本文以水平位移向下游為正，豎向位移向上為正。對應(yīng)圖9中三種情況計(jì)算的壩頂上游側(cè)(點(diǎn)P)水平向位移分別為5.753、5.740和5.859 mm，豎向位移分別為-1.550、-1.553和-1.532 mm。可以看出，是否設(shè)置縫面接觸彈簧對壩頂位移的影響較大，而對于設(shè)置縫面接觸彈簧的情況，粗網(wǎng)格和細(xì)網(wǎng)格的計(jì)算結(jié)果接近，說明網(wǎng)格密度對壩頂位移影響不大。由于在靜力作用下大壩下游預(yù)設(shè)裂縫處于受壓狀態(tài)，變形圖也顯示未施加接觸彈簧的裂縫上下面發(fā)生了嵌入，而施加接觸彈簧的裂縫面處于壓緊狀態(tài)，說明考慮裂縫面的接觸條件是必要的。

圖8 SBFEM靜力計(jì)算模型

圖9 重力壩靜力作用變形圖(變形放大倍數(shù)：600)

5.2 地震斷裂分析結(jié)果仍采用圖8中粗網(wǎng)格和細(xì)網(wǎng)格模型模擬Koyna重力壩地震斷裂過程。地震工況考慮自重+靜水壓力+地震荷載+動(dòng)水壓力，并研究網(wǎng)格密度和裂縫擴(kuò)展步長對大壩地震響應(yīng)和裂縫擴(kuò)展過程的影響。

5.2.1 網(wǎng)格密度對計(jì)算結(jié)果影響 圖10為采用SBFEM粗網(wǎng)格和細(xì)網(wǎng)格模型進(jìn)行非線性斷裂分析得到的壩頂觀測點(diǎn)P的位移時(shí)程，并與不考慮裂縫擴(kuò)展的線彈性分析結(jié)果進(jìn)行對比，由結(jié)果可見：前2.06 s，三個(gè)模型計(jì)算的壩頂水平和豎向位移基本吻合，表明裂縫尚未擴(kuò)展，大壩處于線彈性狀態(tài)；2.06 s時(shí)刻預(yù)設(shè)裂縫起裂，粗、細(xì)網(wǎng)格模型計(jì)算得到的壩頂最大水平向位移分別為-10.06和-10.16 cm，較線彈性分析結(jié)果明顯增大，反映出壩體剛度隨裂縫擴(kuò)展而降低。對比兩種網(wǎng)格計(jì)算的位移結(jié)果可知，壩頂水平和豎向位移時(shí)程相似度極高，且峰值發(fā)生的時(shí)刻一致，說明SBFEM粗網(wǎng)格即可得到滿意的計(jì)算結(jié)果。

圖10 不同模型計(jì)算的壩頂觀測點(diǎn)P的位移時(shí)程

圖11繪制了粗網(wǎng)格和細(xì)網(wǎng)格情況下模擬的壩體裂縫擴(kuò)展軌跡，可看出本文模擬結(jié)果與文獻(xiàn)[8]結(jié)果及振動(dòng)臺(tái)實(shí)驗(yàn)結(jié)果[26]吻合，驗(yàn)證了本文模型的正確性。同時(shí)，兩種網(wǎng)格模擬的裂縫擴(kuò)展軌跡非常接近，綜合圖10中壩頂位移響應(yīng)結(jié)果，說明模型網(wǎng)格密度對動(dòng)態(tài)斷裂分析結(jié)果影響很小，粗網(wǎng)格在動(dòng)力分析中具有足夠的計(jì)算精度。圖12和圖13分別繪制了壩體裂縫擴(kuò)展過程中不同時(shí)刻的變形圖和位移云圖，可以看出，預(yù)設(shè)裂縫起裂后以斜向下的角度迅速向上游擴(kuò)展，在此過程中觀測到地震的往復(fù)作用導(dǎo)致裂縫出現(xiàn)張開和閉合的現(xiàn)象。

圖11 不同網(wǎng)格模型的裂縫擴(kuò)展路徑對比

圖12 Koyna重力壩粗網(wǎng)格模型不同時(shí)刻變形圖(變形放大倍數(shù)：200)

圖13 Koyna重力壩粗網(wǎng)格模型不同時(shí)刻位移云圖(變形放大倍數(shù)：200)

5.2.2 裂縫擴(kuò)展步長對計(jì)算結(jié)果的影響 實(shí)際混凝土起裂后裂縫擴(kuò)展步長受多種因素影響，如加載條件和骨料級(jí)配等，因而有必要探究模型中裂縫擴(kuò)展步長的取值對于斷裂模擬的影響。本次斷裂分析以粗網(wǎng)格模型為例，

分別取0.4、0.8和1.2 m等三種不同的裂縫擴(kuò)展步長，模擬Koyna重力壩的地震斷裂過程。圖14繪制了三種裂縫向上游壩面擴(kuò)展形成貫通裂縫前的變形圖，圖15將不同裂縫擴(kuò)展步長計(jì)算的最終裂縫擴(kuò)展軌跡與文獻(xiàn)[8，26]結(jié)果進(jìn)行了對比。結(jié)果表明，三種裂縫擴(kuò)展步長條件下模擬的大壩裂縫擴(kuò)展軌跡接近；從裂縫擴(kuò)展路徑來看，步長取值越小會(huì)使得裂縫周邊網(wǎng)格更密，裂縫擴(kuò)展軌跡更為精細(xì)。

圖14 不同裂縫擴(kuò)展步長情況壩體變形圖(變形放大倍數(shù)：200)

圖15 不同裂縫擴(kuò)展步長模擬的Koyna大壩裂縫擴(kuò)展軌跡

圖16為采用三種裂縫擴(kuò)展步長計(jì)算的壩頂觀測點(diǎn)P的水平和豎向位移時(shí)程，可看出裂縫擴(kuò)展之前三種情況的壩頂位移時(shí)程吻合，裂縫擴(kuò)展后隨著步長的增大，裂縫張合次數(shù)減少，裂縫向上游擴(kuò)展貫穿壩體用時(shí)越短，對應(yīng)0.4、0.8和1.2 m擴(kuò)展步長取值情況下形成貫穿裂縫用時(shí)分別為4.08、3.36和2.82 s。從圖中也可看出，壩頂水平位移峰值隨著裂縫擴(kuò)展步長的增大而減小，說明裂縫擴(kuò)展步長的取值會(huì)影響壩體的地震響應(yīng)分析結(jié)果。進(jìn)一步對比裂縫張開-閉合一個(gè)周期內(nèi)(如2.3～2.8 s)的位移響應(yīng)，壩頂水平、豎直位移峰值分別為-3.29和-1.25 cm，均來自裂縫擴(kuò)展步長為1.2 m的情況。這是由于擴(kuò)展步長取值越大，壩體剛度降低越快，相同地震作用下對應(yīng)的動(dòng)態(tài)位移響應(yīng)也越大。需要說明的是，混凝土壩動(dòng)態(tài)斷裂過程非常復(fù)雜，每個(gè)時(shí)間步的裂縫擴(kuò)展步長沒有實(shí)測資料可供參考，故本文僅假定了不同的裂縫擴(kuò)展步長以研究其影響，可能與大壩實(shí)際地震斷裂過程存在偏差。

6 結(jié)論

將比例邊界有限元法(SBFEM)和多邊形局部網(wǎng)格重剖分技術(shù)結(jié)合提出了一種全自動(dòng)的重力壩動(dòng)態(tài)斷裂分析模型。該方法繼承了SBFEM降維和半解析描述裂縫尖端奇異應(yīng)力場等優(yōu)點(diǎn)，拓展了SBFEM在斷裂分析中的應(yīng)用，可高精度地求解裂縫動(dòng)態(tài)應(yīng)力強(qiáng)度因子，實(shí)現(xiàn)了重力壩地震斷裂的高效模擬。為追蹤裂縫路徑，僅需增加較少的節(jié)點(diǎn)即可實(shí)現(xiàn)網(wǎng)格重剖分，無須人工干預(yù)，因而具有相當(dāng)高的分析效率。

以Koyna重力壩為例模擬了重力壩的地震斷裂過程，并對模型的網(wǎng)格密度和裂縫擴(kuò)展步長進(jìn)行了敏感性分析。結(jié)果表明：(1)網(wǎng)格密度對壩體地震響應(yīng)和裂縫擴(kuò)展形態(tài)影響很小，且較粗的多邊形網(wǎng)格即可獲得滿意的計(jì)算結(jié)果；(2)裂縫擴(kuò)展步長取值對大壩最終破壞形態(tài)影響不大，但是步長越大裂縫擴(kuò)展越快，從而導(dǎo)致大壩非線性地震響應(yīng)存在差別。綜上，本文方法可成為混凝土壩地震斷裂分析中一種很有競爭力的數(shù)值分析技術(shù)，后續(xù)有望拓展至混凝土界面裂縫以及三維裂縫的動(dòng)態(tài)擴(kuò)展模擬。