混凝土拱壩地震非線性響應中的應變率效應研究

董 雨,王國盛,王進廷,杜修力

(1.北京工業大學 巖土與城市地下研究所,北京 100124;2.清華大學 土木水利學院,北京 100084)

拱壩是一種應用廣泛的壩型,其壩型結構具有優良的力學性能,易解決樞紐泄洪,且容易施工導流。與同級別重力壩相比較,能夠顯著節省筑壩材料,提升經濟效益[1-3]。混凝土作為拱壩主要的建筑材料,是一種典型的率相關材料,其強度受應變率影響十分顯著。對于混凝土應變率效應的研究,最早可以追溯到1917年Abram[4]所進行的混凝土動態單軸抗壓強度的研究,隨后的一百多年里,學者進行了大量關于混凝土率相關強度的研究[5-19]。然而,當前在大壩抗震設計中,通常不考慮混凝土材料的率效應,或者只是在靜強度的基礎上考慮強度的整體放大。例如,應用比較廣泛的是Rapheal[20]利用試驗給出混凝土大壩在地震作用下動強度的建議值,其中混凝土動態抗壓強度較靜強度提高31%,動態直接拉伸強度較靜態直拉強度提高66%,動態劈拉強度提高45%。上述建議值是在確定加載條件下獲得的,即應變速率近似等于5 Hz振動產生的加載速率,但目前已被不分情況地普遍推廣應用于大壩的設計[21]。如2015年NB 35047—2015《水電工程水工建筑物抗震設計規范》[22]中規定,混凝土動態強度標準值可較其靜態標準值提高20%,動態彈性模量標準值可較其靜態標準值可提高50%,動態抗拉強度標準值取為動態抗壓強度標準值的10%。

在實際地震作用下混凝土拱壩的反應中,地震強度、壩體尺度和形式對混凝土的率相關特性均存在較大影響,并且在不同時刻、壩體不同部位混凝土的率效應差別也很大,尤其存在非線性和率效應的耦合影響。混凝土的動態強度若在靜強度基礎上提高20%,無法考慮大壩中混凝土強度在不同時刻和不同部位的變化,原因是沒有采用能夠反映率效應的混凝土本構模型。針對地震循環作用下混凝土材料的損傷特性,Lee等[23]建立了混凝土材料的塑性損傷本構模型,利用該模型開展了Koyna大壩的地震響應分析。陳健云等[24]通過在損傷張量中考慮應變率效應,建立了混凝土材料率相關損傷模型,并開展了300 m級高拱壩地震響應分析,獲得了拱壩混凝土應變率和拉壓損傷的分布規律。肖詩云等[25]將混凝土靜態(Hsieh-Ting-Chen,HTC)模型發展為一致黏塑性模型,對某高拱壩進行了地震反應分析,結果表明應變率對拱壩的應力和塑性應變均有較大影響。其他學者[26-28]相繼發展了混凝土材料率相關本構模型,并開展了高拱壩的地震反應分析。上述研究表明,由于地震荷載在壩體中產生的應變率分布不同,壩體各部位混凝土動態力學性能也會不同。在采用不同的混凝土損傷本構模型時,壩體的損傷分布有所區別,同時大壩的破壞模式與選取的地震波密切相關。在前人所做的數值模擬中,通過采用率相關本構模型開展研究,對考慮應變率和不考慮率影響下壩體特征點的加速度、位移值、壩體拉應力分布和損傷進行了對比分析。從應變率效應產生的本質上來看,地震荷載在壩體中產生的應變率隨時間和空間不斷變化,變化的應變率會引起混凝土動態強度變化,最終導致壩體的加速度、位移、應力分布情況、損傷值等是連續變化的過程,而非簡單的一個狀態。

本文采用有效硬化函數和動態S型強度準則對混凝土塑性損傷模型進行改進,利用S-CDP模型開展了大崗山拱壩有限元模擬。從壩體應變率的角度進行研究,選取壩體上具有代表性的一些單元,將特征單元的應變以及應變率隨地震時間的變化曲線提取出來,通過直接對不同特征單元的應變率曲線進行分析,得到壩體不同部位應變率隨時間的變化過程。再從宏觀層面上對由于動強度的不同引起的損傷情況進行研究,分析地震作用下混凝土拱壩的應變率分布、相應的動態增長因子(dynamic increase factors,DIF)的變化情況以及應變率效應帶來的影響。

1 混凝土損傷塑性模型的改進

ABAQUS有限元軟件中內嵌的混凝土損傷塑性(concrete damage plasticity,CDP)模型[29]常被作為混凝土的材料模型使用,原始CDP模型需要用戶自定義動態強度隨應變率變化規律,然后將單軸拉壓應力-應變曲線作為非彈性應變率的函數來考慮應變率效應。動強度定義不合理很容易導致數值計算不收斂。針對這一問題,本文對原始CDP模型進行了改進。首先采用在有效應力空間中確定的多軸硬化函數代替名義應力空間中的單軸硬化函數,保證三維塑性變形的合理計算;其次選用連續光滑、能反映真實動強度的動態強度S準則[30]考慮混凝土的應變率效應,在子程序中使用S準則來計算應變率相關應力-應變關系,提高計算的收斂性,同時避免慣性效應的多次考慮,從而減小計算結果與實際結果的偏差。具體的改進方法闡述如下。

原始CDP模型為雙標量塑性損傷本構模型,一個損傷變量用于拉伸損傷,一個損傷變量用于壓縮損傷。通過損傷變量建立了名義應力與有效應力的關系

(1)

(2)

式中：E0為初始無損傷彈性剛度張量;ε為總應變張量;εp為塑性應變張量。

式(1)中ω為包含壓縮損傷和拉伸損傷的損傷變量

ω=1-(1-ωt)(1-ωc)

(3)

式中：ωt為拉伸損傷變量;ωc為壓縮損傷變量。

有效應力表示的屈服函數是得到塑性演化的關鍵,CDP模型的屈服函數表達式如下

cc(κ)

(4)

(5)

在屈服函數式(4)中,cc(κ)和β(κ)為兩個硬化函數。其中cc(κ)為壓縮硬化,β(κ)描述壓縮和拉伸的協同硬化,其表達式為

(6)

式中：ct(κ)為拉伸硬化函數;κ為硬化參數。

1.1 有效硬化函數

(7)

采用多軸應力條件下的統一硬化/軟化函數[33]來描述名義應力空間中的硬化函數,它能夠考慮應力狀態對硬化/軟化行為的影響,其方程如下

(8)

式中：A和D為形狀參數;xi=εpd/εpd0為硬化參數;εpd為等效塑性剪切應變;εpd0為峰值應力下εpd的值。εpd的表達方式為

(9)

式中,εpi為塑性主應變,本文規定拉為正,壓為負。

將式(8)代入式(7),得

(10)

式(10)為多軸條件下的有效硬化函數,對于單軸拉伸條件,式(9)中的εp2=εp3=0,εp1>0;對于單軸壓縮條件,式(9)中的εp1=εp2=0,εp3<0。因此,將上述應變條件代入式(10),多軸硬化規律可退化為單軸拉伸和單軸壓縮硬化規律,進而可以確定屈服函數式(4)中的cc(κ)與β(κ)。當CDP模型采用由有效應力空間中確定的硬化規律時,S-CDP模型能夠更加合理地描述混凝土材料在有效應力空間中的屈服面演化規律。

1.2 應變率相關規律

對于混凝土材料的動態單軸強度準則,最具代表性的準則為CEB-FIP規范建議的雙線性單軸動態強度準則[34],該準則的方程是分段函數,對于分段表達形式的動態單軸強度準則,準則中存在拐點,拐點位置的確定缺乏理論依據。當前大多數既有的動態強度準則都是發散的,無法反映混凝土的極限動強度,并且也沒有區分真實動強度與宏觀動態承載力,導致計算結果過高估計混凝土材料的真實動強度[35-36]。在本文將引入非線性動態單軸S準則,該準則是基于混凝土材料的動態物理機制建立的,能夠合理反映混凝土材料的真實動強度和試驗表現出的極限動強度,S準則的材料參數較少且具有明確的物理意義,其強度曲線連續光滑,不存在明顯的拐點。S準則的方程如下

(11)

2 大崗山拱壩數值模型及地震作用下的損傷規律

2.1 大崗山拱壩模型參數與建模過程

大崗山拱壩壩型為混凝土雙曲拱壩,最大壩高為210 m,壩頂高程為1 130 m[38]。大崗山壩體-地基有限元模型如圖1所示。數值模型設置了與實際工程相同的28條壩體橫縫[39],在進行分析的過程中,考慮了鍵槽咬合作用的高拱壩非線性行為[40]。拱壩橫縫在地震荷載作用下張開、閉合的非線性力學行為引入了動接觸邊界以及限制接觸面切向滑移的彈簧單元[41]。并采用黏彈性人工邊界模擬無限地基輻射阻尼效應[42]。

(a) 壩體-地基

參照同類工程經驗,本文的計算模型選取的地基范圍為：下游面地基2倍壩高,上游面地基、左右岸及自壩底向下均為1.5倍壩高。有限元模型中,采用S-CDP模型作為混凝土的材料模型;壩體混凝土單元類型采用八節點六面體線性縮減積分實體單元(C3D8R),模型總的單元數量為121 122個,其中壩體單元數為79 638個,地基單元為41 484個,節點數為148 551。

分別采用“Static”分析步計算初始靜態荷載的作用,采用“Dynamic/Implicit”分析步計算地震荷載的作用,初始增量步時長設置為0.01 s。壩體混凝土的材料參數為：混凝土密度2 400 kg/m3,彈性模量較其靜態標準值提高50%取為36.65 GPa,泊松比0.17,膨脹角45°,α=0.14,形狀參數A=2.5,D=2.0。S準則中強度參數的取值：動態拉伸條件下ξ=1,Fmax=10.2,u0=-1,F0=0.5;動態壓縮條件下ξ=1.2,Fmax=3.72,u0=1.95,F0=0.5。基巖假定由兩部分構成：遠離壩體的線彈性基巖與靠近壩體的非線性基巖。線彈性基巖與非線性基巖相同的材料參數為：基巖密度2 650kg/m3,彈性模量20 GPa,泊松比0.25。非線彈性基巖采用Drucker Prager Harding模型,摩擦角65.1°,流變應力比取1,膨脹角取65.1°。作用荷載包括壩體自質量、上游庫水靜水壓力、溫度荷載、上游庫水動水壓力和地震荷載。其中,動水壓力采用Westergaard附加質量法考慮[43]。

大崗山拱壩所處地區地震基本烈度為Ⅷ度,設防標準按100 a設計基準期超越概率為2%(重現期5 000 a)的基巖水平加速度,其設計加速度水平向為0.557 5g,豎向加速度取水平方向的2/3[44]。為突顯應變率對拱壩地震非線性響應的影響,本文考慮地震動的不確定性[45],從實測地震數據庫中選取3組地震動時程曲線：第一組為等比例放大1.271倍得到峰值地面加速度(peak ground acceleration,PGA) 0.557g的地震動時程曲線。第二組等比例放大1.821倍得到PGA=0.663g地震動時程曲線;第三組為等比例放大1.025倍得到PGA=0.836g的地震動時程曲線;XYZ3個方向的地震波加速度時程曲線如圖2所示。

(a) PGA=0.557g的地震動時程曲線

2.2 3組地震作用下壩體的損傷規律

利用S-CDP模型作為混凝土的材料模型,分別計算了在第一組、第二組、第三組地震動荷載作用下壩體的動力響應,其損傷分布如圖3～圖5所示。由圖3～圖5可知,壩基面的損傷比較嚴重,壩體下游面出現了大面積的損傷,損傷結果與其他學者[46]論文中的強震下高拱壩的損傷開裂結果具有一定的相似性,說明了本文計算結果的合理性。在第一組地震的損傷結果中,上游面靠近壩頂處也出現了小范圍的損傷。兩組地震作用下壩體損傷最嚴重的部位均出現在壩基,原因是壩基部位的應力集中所導致。

(a)下游面

(a)下游面

(a)下游面

在3組地震作用下得到的絕對值最大的主應變如圖6所示,將以符號εmpa來表示絕對值最大的主應變。選取7個特征單元,編號為：L1,L2,R1,R2,F1,F2,M1。

(a) 第一組地震動下的分布情況

壩體的εmpa主要分布在壩踵部位及壩體下游面的中部,結合損傷分布的情況來看,壩體在承受地震荷載作用時,最容易出現損傷的部位為壩基以及壩體下游面的中部區域。

第一組地震動、第二組地震動和第三組地震動荷載作用下壩體7個特征單元的εmpa時程曲線如圖7所示。圖7表明,選擇的特征單元的拉應變要大于壓應變,而壩體混凝土的抗拉強度低,因此出現的損傷都是拉伸損傷,沒有壓縮損傷。

(a) 第一組地震動

對圖7的εmpa的時程曲線進行分析,第一組地震動下壩基處兩個特征單元F1和F2的εmpa最大分別達到了0.004 5和0.006 5。除壩基外壩體的其他特征單元的εmpa遠小于特征單元F1和F2的值,都分布在0.000 5以下。第二組地震動荷載作用下,壩基處的兩個特征單元F1和F2的εmpa最大分別達到了0.007和0.007 5。除壩基外壩體的其它特征單元的εmpa遠小于F1和F2單元的值,都分布在0.002以下。在第三組地震動荷載作用下,壩基處的兩個特征單元F1和F2的εmpa最大分別達到了0.002 2和0.006,除壩基外壩體的其他特征單元的εmpa遠小于特征單元F1和F2的值,都分布在0.001以下。

3 混凝土拱壩的應變率效應分析

為研究在地震荷載作用下混凝土拱壩的應變率效應,采用S-CDP模型,對大崗山混凝土拱壩進行了第一組地震動、第二組地震動與第三組地震動下的非線性地震響應分析,得到其應變率變化規律和相應的DIF變化規律。

3.1 第一組地震動作用下的應變率規律

基于DIF的定義(動強度與靜強度的比值),通過應變率時程曲線獲得第一組地震動作用下特征單元的DIF時程曲線,如圖8所示。

(a) 拉伸動態增長因子

在第一組地震動作用的情況下,拉伸DIF變化范圍為1.00～1.21,壓縮DIF變化范圍為1.000 0～1.003 3。混凝土的應變率效應增大了壩體動強度,以拉應變產生的動強度增長為主,DIF的數值分布表明動強度增長范圍大致分布在0%～21%。

3.2 第二組地震動作用下的應變率規律

對圖7中εmpa的時程曲線微分,得到特征單元拉應變率和壓應變率隨時間的變化規律如圖9所示。圖9(a)表明：在整個地震過程中拉應變率幅值最小的是R1單元的曲線,該單元位于壩體的中部偏右岸位置;拉應變率幅值最大的是F2單元的曲線,該單元位于壩體壩基部位。拉應變率最大峰值出現的時間段分布在地震荷載作用的5～7 s,該時間段也是地震加速度時程曲線最大峰值出現的區間。圖9(b)表明：壓應變率幅值最小的是R2單元的曲線,該單元位于壩體的中間部位;壓應變率幅值最大的為F2單元的曲線。

(a) 拉伸應變率

分析圖9中應變率的變化規律,壩體出現拉壓應變率幅值最大的位置與損傷最大的位置是一致的。然而,拉應變率與壓應變率幅值最小的單元則出現在壩體的不同部位。此外,壩體不同位置處單元應變率的幅值和出現幅值的時刻都是不同的。通過采用S-CDP模型,計算得到壩體多處單元的拉壓應變率曲線,直觀地展示了壩體不同部位混凝土應變率在時間和空間上的連續變化過程。

基于DIF的定義(動強度與靜強度的比值),計算得到所選單元的DIF隨時間的變化規律,如圖10所示。

(a) 拉伸動態增長因子

在第二組地震動作用的情況下,壩體的拉應變率變化范圍為0～1.6×10-2s-1,最大拉應變率只在6.2 s時出現一次,其余時刻壩體混凝土的拉應變率基本上都位于0.8×10-2s-1以下。壩體的壓應變率范圍為0～0.25×10-2s-1,其余時刻壩體混凝土的壓應變率基本上都位于0.15×10-2s-1以下。在圖10(a)中,拉伸DIF幅值最小的是R2單元的曲線,拉伸DIF幅值最大的是F2單元的曲線。在圖10(b)中,壓縮DIF幅值最小的是R2單元的曲線,壓縮DIF幅值最大的是F2單元的曲線。拉伸DIF變化范圍為1.00～1.23。壓縮DIF變化范圍為1.000 0～1.002 5。拉伸DIF的最大值遠大于壓縮DIF的最大值,但是應指出拉伸DIF最大值1.23只在6.2 s時出現了一次,出現的部位為壩基處的F2單元。

綜上數據,采用S-CDP模型能夠很好地考慮不同應變率下混凝土動強度的變化,且不同應變率下的動強度和已有的工程經驗是一致的,證明了模型的合理性。在第二組地震動荷載作用下,混凝土的應變率效應對壩體動強度增長起到了一定的作用,以拉應變產生的動強度增長為主,DIF的數值分布表明動強度增長范圍大致分布在0%～23%。

3.3 第三組地震動作用下的應變率規律

采用與3.2節相同的方法,通過應變率時程曲線獲得第三組地震動作用下特征單元的DIF時程曲線,如圖11所示。

在第三組地震動作用的情況下,壩體的拉應變率變化范圍為0～1.1×10-2s-1,壩體的壓應變率范圍為0～0.27×10-2s-1,在整個地震動荷載作用的時間內,壩體混凝土的拉應變率基本上都位于0.4×10-2s-1以下,壩體混凝土的壓應變率基本上都位于0.2×10-2s-1以下。在圖11(a)中,拉伸DIF幅值最小的是R1單元的曲線,拉伸DIF幅值最大的是F2單元的曲線。在圖11(b)中,壓縮DIF幅值最小的是M1單元的曲線,壓縮DIF幅值最大的是F2單元的曲線。拉伸DIF變化范圍為1.00～1.19,壓縮DIF變化范圍為1.000 0～1.002 5。拉伸DIF的最大值遠大于壓縮DIF的最大值,但是最大的拉伸DIF值1.19只在5.4 s時出現了一次,出現的部位為位于壩基處的F2單元。

綜上數據,在第三組地震動作用下,混凝土的應變率效應增大了壩體動強度,以拉應變產生的動強度增長為主,DIF的數值分布表明動強度增長范圍大致分布在0%～19%。

4 率相關強度對壩體地震響應的影響

為了分析率相關強度對壩體地震響應的影響,本章考慮在第二組和第三組地震動荷載作用,分析3種工況：工況I,采用率無關CDP模型,按混凝土材料靜態強度,彈性模量較其靜態標準值提高50%計算的拱壩非線性響應;工況II,按照現行的NB 35047—2015《水電工程水工建筑物抗震設計規范》,直接將壩體混凝土材料靜態強度提高20%,動態彈性模量標準值較其靜態標準值提高50%,采用率無關CDP模型計算拱壩非線性響應;工況III,采用S-CDP模型作為混凝土材料模型,動態彈性模量標準值較其靜態標準值提高50%,計算拱壩的動力非線性響應。3種不同的工況,均同時考慮壩體自體質量、靜水壓力、動水壓力、溫度荷載和地震荷載作用,第二組地震動荷載壩體損傷情況如圖12所示,第三組地震動荷載壩體損傷情況如圖13所示。

(a) 工況Ⅰ

(a) 工況Ⅰ

圖12和圖13中分別展示了第二組和第三組地震動作用下的3種工況壩體下游面的損傷情況。從圖12、圖13可知,兩組地震動荷載的3種工況壩體下游面都出現了損傷情況但均未形成自下游壩面至上游壩面的貫穿損傷。對比工況I、工況II、工況III的損傷響應,可以發現,在同一地震動荷載作用下,工況Ⅰ和工況Ⅲ壩體的損傷區域基本一致,都出現在壩基以及壩體的中上部位,但是工況Ⅲ中的壩體損傷程度比工況Ⅰ的損傷程度小,這說明在地震荷載的作用下,壩體混凝土采用S-CDP模型計算得到的損傷比采用率無關CDP模型計算得到的損傷會減小,率相關強度對最終的損傷結果產生了影響。

對圖12(b)和圖12(c)的損傷云圖及損傷因子進行對比,對圖13(b)和圖13(c)的損傷云圖及損傷因子進行對比,可以看出工況Ⅱ的壩體損傷區域與損傷因子較工況Ⅲ整體偏小,但是在實際工程中混凝土拱壩的地震響應存在應變率效應,其壩體的損傷區域及損傷因子與在靜強度的基礎上考慮強度的整體放大20%的情況相比時會更大,這說明采用S-CDP模型時混凝土的動強度增長無法一直達到20%的程度,也說明在進行水工建筑物抗震計算時,若將整個壩體的混凝土的動強度提高20%,可能會造成損傷結果與實際情況有所偏差。

5 結 論

通過引入有效硬化函數和S型率相關強度準則,考慮混凝土材料的三維塑性變形行為和應變率效應,對混凝土塑性損傷本構模型進行了改進。使用S-CDP模型,對大崗山混凝土拱壩進行了非線性地震響應分析。分別計算了PGA=0.557g,PGA=0.663g,PGA=0.836g的地震動荷載作用下拱壩的動力響應,分析了壩體的應變分布規律、應變率的分布規律及相應的DIF變化規律、證明了S-CDP模型的合理性。主要結論如下：

(1) 在地震荷載作用下,拱壩混凝土的壓應變率效應不明顯;然而,壩基和部分壩體區域產生了較大的塑性變形,拱壩混凝土的拉應變率效應明顯。此外,通過分析壩體中特征單元的應變以及應變率隨時間的動態變化規律,發現采用S-CDP模型建立的數值模型能夠合理反映混凝土的應變率效應。

(2) 通過在本構模型中考慮混凝土率效應,計算得到的大崗山拱壩壩體的拉應變率最大峰值為1.6×10-2s-1,相對應的拉伸DIF為1.23,僅在壩基處的某瞬間出現了一次,除此外拉應變率變化較大,并不能時刻保證壩體混凝土的動強度提高20%及以上,地震荷載作用下壩體混凝土的壓應變率效應不顯著。在水工建筑物抗震設計時,若將混凝土的動強度提高固定的20%,可能會造成計算的損傷結果與實際的損傷有所偏差。動強度提高固定的20%計算得到的損傷結果相較于本文中采用S-CDP模型計算得到的損傷結果偏大。因此在水工建筑物抗震設計時,推薦選用合適的率相關本構模型進行計算,以盡可能的保證計算的精確度。