核電廠考慮結構參數不確定的地震易損性分析

秦苗珺, 趙衍剛,2, 盧朝輝

(1. 北京工業大學 城市建設學部, 北京 100124; 2. 神奈川大學 工學部建筑學科, 日本 橫濱 221-0806)

0 引言

核電廠發生事故后,對社會和環境將產生難以估量的影響,因此,核電廠的安全運營是關乎經濟發展與社會安定的根本。核電站設計建造過程中應考慮各種不確定因素對其安全性的影響,并針對這些因素采取相應有效的措施。自從福島核事故以來,越來越多的研究人員對既有核電站和新建核電站重新進行安全評估。由于同一震級下不同結構的地震烈度有明顯區別,不同的烈度主要取決于建筑的結構材料,施工質量等因素,其中結構的材料屬性可反映施工質量的優劣,影響著結構的整體性能。因此,材料參數不確定性對結構抗震性能的影響仍不可忽略。

在核電廠的抗震分析中,對核電站在超過基準地震動作用下的有效評估是保證核電安全的關鍵。因此,核電站的地震概率風險評估(Seismic Probability Risk Assessment,SPRA)已廣泛應用于新建或已建的核電站系統的安全評估中[1]。而地震易損性分析是SPRA中至關重要的部分,是評估在地震作用下結構可靠性的關鍵步驟,為抗震設計與風險評估提供指導依據。因此,準確獲得結構的易損性是地震風險評估的前提條件。目前,美國所提出的多種地震易損性分析方法已被世界許多國家所采用[2],主要為以下三種方法:Zion法[3]、地震安全裕度法(Seismic Safety Margin Research Program,SSMRP)[4]和BNL(Brookhaven National Laboratory)[5]方法。其中Zion法與SSMRP法都通過專家評估與經驗判斷確定,結果具有較大的主觀性。在BNL方法中,通過隨機振動理論與極值理論獲得結構的最大反應分布,而參數不確定則是通過拉丁超立方抽樣進行考慮。Whittaker等[6]綜述了隔震核電廠地震易損性的發展。Zhao 等[7]對屏蔽廠房進行了地震易損性評估,考慮了流固耦合對其影響,同時采用三種不同方法建立了核電廠系統的易損性曲線。

在實際工程結構中除地震動激勵的不確定外,還存在外部環境、材料屬性和人工干預等不確定因素,而多數核電廠的風險評估中未考慮參數不確定的影響[8]。研究表明,參數不確定對結構可靠性分析有不利的影響[9]。早期是通過蒙特卡洛模擬(Monte Carlos Simulation,MCS)[10]考慮參數不確定,但由于模擬方法對于小失效概率事件計算量龐大,并且獲得結構的響應函數計算時間成本較大,難以直接使用。隨后許多學者針對該方法進行了改進,發展了不同的抽樣方法,N Yun等[11]提出采用一組模型的輸入輸出樣本估計整體可靠性靈敏度指標的方法,并結合子集模擬以減小計算量;Alvarez等[12]采用隨機集理論計算失效概率的上下界,其中變量之間的相關性采用Copula函數表示,該方法對于低維空間有較好的適用性。Alban等[13]提出了一種高效處理高維、小失效概率問題的子集模擬方法,選擇合理的中間失效事件,將較小的失效概率表達為一系列較大失效概率的乘積,并利用馬爾科夫鏈模擬(Markov Chain Monte Carlo,MCMC)方法生成條件樣本計算失效概率,該方法對變量維數、極限方程形式等均沒有限制,適用于非線性較高的小失效概率可靠性問題,但對于中間失效事件的選取及生成樣本的相關性對計算結果有較大的影響,為減小其影響,發展了自適應重要抽樣法等[14-15]。其次以展開法為主的矩法等是考慮參數不確定的有效方法[16]。Zhao Y G等[17]提出高階矩方法,對已有方法進行了簡化,同時具有較高的精度。一次二階矩法(First-Order Second-Moment,FOSM)與MCS方法相比因計算簡單已被廣泛應用于各種結構,如混合結構、鋼結構、橋梁結構的易損性分析中[18-19]。蔣亦龐等[20]基于FOSM方法分析了無筋砌體結構在參數不確定下的地震易損性,結果表明,參數不確定對易損性的影響不可忽略。

本文以混凝土的材性參數為隨機變量,采用ANSYS有限元建模,通過與試驗結果對比驗證,表明有限元模型可作為易損性分析中的模擬工具。在此基礎上結合增量動力法(Incremental Dynamic Analysis,IDA)[21]與FOSM考慮核電廠在參數不確定下的地震易損性,得到核電結構的地震易損性曲線。結果表明,參數不確定對核電結構的地震易損性具有顯著的影響。

1 核電廠結構模型與驗證

1.1 核電站的有限元模型

本文利用有限元軟件ANSYS的交互仿真平臺,對CAP1400核島結構及周邊附屬廠房建模。其尺寸如圖1所示。本文主要考查核島結構在地震作用下結構關鍵點的動力響應,因此在建模時對實際的復雜結構進行適當的簡化,不考慮屏蔽廠房內部鋼制安全殼的影響。因結構的厚度方向尺寸相對較小,屏蔽廠房和輔助廠房采用殼單元,選用SHELL181薄殼單元進行建模,基礎底板采用實體SOLID65單元。保證混凝土材料全曲線有下降段,本文采用了多線性隨動強化(MKIN)模型考慮混凝土進入塑性的性能(圖2)。為滿足規范要求,采用C50混凝土,抗壓強度fc=23.1 MPa,密度rc=2 500 kg/m3,彈性模量Ec=3.45×104MPa,泊松比uc=0.18。混凝土的應力應變曲線采用Kent-Park模型[22]。屏蔽廠房和輔助廠房之間、結構和土體之間采用MPC型綁定約束。場地巖土力學參數如下:rs=2 600 kg/m3,Es=10 030 MPa,剪切波速vs=1 250 m/s,泊松比us=0.32,內聚力cs=0.9 MPa,摩擦角fs=42.8°,選用ANSYS中自帶的Drucker-Prager模型。土體周圍邊界參考劉晶波等[23]和Deeks等[24-25]提出的黏彈性人工邊界,如圖3所示。其中,彈簧阻尼器單元的參數按式(1)和式(2)確定。

圖1 核電站結構模型平面圖(單位:mm)Fig.1 Plan of the nuclear power plant structure model (Unit:mm)

圖2 混凝土線性隨動強化模型Fig.2 Linear kinematic hardening model of concrete

圖3 三維黏彈性人工邊界示意圖Fig.3 Diagram of the 3D viscous-elastic artificial boundary

法向:

(1)

切向:

(2)

根據以上結構參數建立的核島結構有限元模型如圖4所示。

圖4 核電站的整體模型Fig.4 Overall model of the nuclear power plant

1.2 模型驗證及對比

為驗證有限元模型的準確性,建立與振動臺試驗相一致的模型進行對比。試驗模型采用1∶40縮尺比例,縮尺后的模型尺寸為2 285 mm×1 058 mm×1 188mm(長×寬×高),屏蔽廠房尺寸1 200 mm×2 194mm(直徑×高),底板厚160 mm。附屬廠房和屏蔽廠房部分模型如圖5所示。地震動輸入選用安縣地震記錄。安縣地震動是2008年汶川地震中觀測到的地震加速度記錄[26],其中水平X方向的峰值加速度為2.99 m/s2,持續時間為50 s,時間步長為0.008 s。地震動加速度時程如圖6所示,相應的加速度反應譜如圖7所示。峰值加速度采用0.3g。選取屏蔽廠房沿標高從上到下7個關鍵點進行對比。圖8展示了安縣地震動下個測點峰值加速度的試驗結果與有限元分析結果的對比。

圖5 核電站試驗的縮尺模型與振動臺試驗的 整體模型Fig.5 Scale model of the nuclear power plant and overall model of shaking table test

圖6 輸入的地震動及加速度時程Fig.6 Input ground motions and acceleration time history

圖8 安縣地震動下各測點峰值加速度的試驗結果與有限元分析結果的對比Fig.8 Comparison between test result and finite element analysis result of peak acceleration of each measuring point under Anxian ground motion

從圖8中可看出,有限元模型在一定程度上反映了結構的真實動力特性。從X方向可看出,振動臺試驗中,核島結構下部加速度響應較有限元結果偏小,數值模型的整體加速度響應從底部到頂部變化范圍小,有限元模型沿高度方向基本呈線性變化,表明結構處于彈性狀態,整體剛度比試驗模型偏大;結構頂部的峰值加速度比有限元模型的偏大,其主要原因可能是X方向的結構剛度沿層高分布不均勻,同時在試驗中土體屬性由于振動過程有所改變,與結構之間的接觸邊界也相應發生變化。而Y、Z向峰值加速度的數值分析結果與試驗相比偏大,可能由于模型中兩方向的整體剛度較大的原因,試驗中土體的整體剛度未達到預期效果,而且由于土層在振動臺試驗過程中整體形態有所變化,與結構之間的相互作用減弱。因此,在軟弱地基下造成結構響應偏小。從上圖表明該有限元模型在一定程度上可等效為實際的核電廠結構。

通過將上述縮尺模型按照相應比例放大為原模型尺度,場地土采用相同的土體參數,人工邊界采用同樣的黏彈性人工邊界,并利用該模型進行易損性分析。

2 地震動參數及結構極限狀態

2.1 地震動與地震強度指標選取

通過已有研究表明,在IDA中選取20條地震記錄足以考慮地震動輸入的不確定性[4,9]。本文依據AP1000反應譜(圖9)從美國太平洋地震研究中心PEER的強震數據庫中選取20條實測地震記錄,震級分布在5.2～7.49,峰值加速度分布在0.08g～0.76g,如表1所列。

表1 20條天然地震動記錄

圖9 AP1000設計反應譜Fig.9 Design response spectrum of AP1000

地震易損性分析中常用的地震動強度參數(Intensity Measure,IM)一般取為結構基本周期對應的加速度譜值Sa(T1,5%)(T1為結構的基本周期)或峰值加速度 (Peak Ground Acceleration,PGA)。由于結構參數不確定會導致結構基本周期成為一個不確定變量,采用Sa(T1,5%)作為IM參數會使分析變得復雜,因此,選取PGA作為IM參數。在進行IDA分析時,分別將20條地震記錄的PGA調整為0.05g～0.8g(間隔為0.05g),對結構進行有限元分析。

2.2 核電結構不確定參數

由于核島結構以混凝土材料為主,因此選用混凝土的密度rc,彈性模量Ec,泊松比uc,抗拉強度ft作為隨機變量。由于目前針對核電結構不確定性參數之間的相關性研究仍不充分,現有文獻難以獲得參數間的相關系數,因此本文采用簡化分析處理,假設各參數之間相互獨立。同時,本文也不考慮模型參數在結構空間分布上的不確定性,假設四個隨機變量的概率信息如表2所列。

表2 結構不確定性參數的概率信息

2.3 結構極限狀態的定義

本文參考了Crowley等[27]提出的基于失效模式反演結構在不同極限狀態的抗震能力的方法,通過基于變形與應力混合控制的方法定義極限狀態,假設結構的最大拉應變達到混凝土的極限拉應變作為結構輕微破壞的極限狀態(Limit State,LS)LS1;中等破壞為LS2,定義為混凝土最大壓應變達到極限壓應變的一半;嚴重破壞定義為混凝土達到峰值壓應力所對應的極限狀態LS3。

文中選取屏蔽廠房結構的頂部位移作為結構抗震性能指標,通過PUSHOVER分析確定核島結構的三個極限狀態對應的頂部最大位移,如表3所列。

表3 頂部最大位移限值

3 考慮結構參數不確定性的地震易損性分析方法

3.1 地震易損性模型

結構地震易損性是指在不同強度的地震作用下,結構超過某一特定極限狀態的失效概率,從宏觀上反映了地震動強度與結構的損傷程度之間的關系,是評估結構抗震性能水平的關鍵組成部分。通常結構的地震易損性計算模型可表達為在某一確定的地震強度IM下,地震需求D達到或超過結構抗震能力C的條件概率,即,

FG(x)=P[D≥C|IM=x]

(3)

式中:FG(x)為易損性函數,通常采用對數正態分布作為易損性概率模型[28-29],即認為地震需求D與抗震能力C均服從對數正態分布,對于某一特定極限狀態LS,地震易損性可表示為:

FG(x)=P[LS|IM=x]=

(4)

式中:Φ[·]是標準正態變量的累積分布函數;Ci為極限狀態LSi下結構抗震能力限值;mG與βG分別為結構地震易損性函數的中位值和對數標準差。當不考慮參數不確定時,βG=βRTR,表征僅有地震動不確定性的影響,可根據16%分位曲線與84%分位曲線,按照式(5)計算得到[30]:

mR,16%=mGexp(-βRTR)
mR,84%=mGexp(βRTR)

(5)

式中:mR,16%表示16%分位曲線;mR,84%表示84%分位曲線。

通過式(5)可得到相應的地震易損性曲線。地震易損性曲線直觀反映了結構的抗震性能,曲線的中心點取決于結構的抗震極限承載能力,而其形狀主要由結構抗震性能的不確定決定,如圖10所示。當結構參數不確定性較大時,不考慮不確定性所帶來的影響容易高估結構的可靠性,使得結構偏于不安全。因此,在參數不確定的條件下,有必要結構參數對結構地震易損性的影響。

圖10 地震易損性曲線的概率特征Fig.10 Probability characteristics of seismic fragility curve

3.2 基于FOSM-IDA考慮參數不確定的易損性分析方法

本文采用FOSM考慮參數不確定的影響。假設結構存在隨機變量X=[X1,X2,…,Xn],諸如結構材料的彈性模量,強度等。易損性函數如式(4)中,由于參數不確定的影響,mG與βG均體現為隨機變量的形式,即考慮參數不確定的易損性模型轉化為對模型參數的估計,以下將通過FOSM確定對數均值mG與對數標準差βG。

假設結構的功能函數Z為隨機變量X的函數Z=G(X),FOSM方法是通過將函數在隨機變量X的均值μX處近似展開為一階泰勒級數形式估計相應的均值和方差[19],如式(6)和式(7)所示:

μZ≈G(μX)

(6)

(7)

(8)

當隨機變量相互獨立時,式(7)可寫為:

(9)

據文獻[31]通過數值模型與式(9)即可得到各隨機變量在各個取值(即平均值加減一倍標準差)下的地震需求與參數不確定下的對數標準差。綜合考慮地震動不確定和結構參數不確定后,此時新定義的易損性函數的總對數標準差βTOT如式(10)所示:

(10)

式中:βRTR與βMOL分別表示只考慮地震動不確定與結構參數不確定的對數標準差,可分別根據式(5)和式(9)計算得到。

4 核電結構易損性分析結果

選用20條實際地震動記錄作為輸入(表1),對核電站結構進行IDA分析,并按照上述方法考慮不確定參數下的地震易損性分析,得到各峰值加速度下不同地震記錄的結構最大響應值,通過采用基于IM準則得到16%、50%和84%分位數曲線如圖11所示。

圖11 16%、50%和84%分位數曲線Fig.11 16%, 50%, 84% quantile curves

由圖11可知,以上三條分位數曲線在起始階段有明顯的線彈性段,表明結構響應與PGA是線性關系,隨著PGA的增加,曲線斜率明顯減小,結構響應出現非線性。

通過對IDA結果進行統計回歸后,得到未考慮結構參數不確定的地震易損性曲線,同時基于FOSM計算結構參數不確定下核電結構不同極限狀態的地震易損性曲線,如圖12所示。其中NON-表示只考慮地震動不確定時的計算結果,CON-表示同時考慮地震動和結構參數不確定性時的計算結果。

圖12 結構易損性曲線的對比Fig.12 Comparison between structural fragility curves

根據地震易損性曲線的特性可知,曲線的傾斜程度反映了不確定性對地震易損性的影響程度。當PGA處于0.7g～1.4g時,達到LS2極限狀態的失效概率在0.2～0.8范圍內,表明結構響應存在較大的離散性,并且在LS3極限狀態表現更加明顯。通過圖12可看出考慮參數不確定與未考慮參數不確定的兩條易損性曲線不重合,考慮參數不確定的易損性曲線斜率小于未考慮參數不確定的易損性曲線,未達到地震強度中位值前,考慮參數不確定后計算失效概率高于未考慮參數不確定的情況,并且隨著破壞等級提高,其相差越大。當在LS2破壞等級,考慮地震動輸入不確定與參數不確定的對數標準差分別為βRTR=0.235和βMOL=0.174(表4),兩者比值βMOL/βRTR為0.738。可見,結構材料參數不確定對核電站廠房結構的地震易損性的影響不能忽略。

表4 結構參數不確定對地震易損性的影響

為進一步研究以上四個隨機變量對核電結構的地震易損性的影響,對不確定參數進行敏感性分析。在上述計算的基礎上,隨機參數分別取均值增減1倍和2倍標準差進行地震易損性分析。此時,每個參數下包括5個計算工況:μXi,μXi±σXi,μXi±2σXi。以易損性曲線的地震強度中位值mR作為評價指標,對每個結構參數下的5個計算工況結果進行歸一化處理:

(11)

結果如圖13所示。從圖中可看出,除彈性模量外,其他參數均在0.87～1.08范圍內變化,變化幅度為0.21。而彈性模量的歸一化結果的變化范圍是0.76～1.23,幅度為0.46,是前者的兩倍之多。表明材料彈性模量的敏感性遠高于其他三個參數,是結構地震易損性分析的控制性參數。

圖13 地震強度中位值隨各結構材料參數的變化規律Fig.13 Variation of the median value of seismic intensity with different structural material parameters

5 結 論

以核電站廠房為研究對象,采用IDA方法計算地震易損性,并基于FOSM考慮結構材料參數的不確定對其影響,主要結論如下:

(1) 在核電廠結構地震易損性分析中,需要同時考慮地震動與結構材料參數兩類不確定的影響,并且隨著破壞等級的提高,結構不確定性參數UI易損性的影響程度也越大。

(2) 結構參數的不確定相比于地震動不確定是不可忽略的,結構參數不確定的標準差占到地震動不確定標準差的60%,因此有必要在核電站易損性分析中同時考慮兩種不確定的影響,同時,在核電站的四種不確定參數中,以彈性模量對地震易損性的影響最為顯著,是參數不確定中的控制因素。

(3) FOSM-IDA便于計算,可作為工程中一種計算地震易損性的有效方法,通過考慮參數不確定的影響,能夠有效解決原有方法在低失效概率情況下高估結構抗震性能的問題,更加準確地評估結構的易損性。