預應力混凝土安全殼內壓易損性及性能評估

金松， 李忠誠， 貢金鑫， 董占發， 藍天云

(1.大連理工大學 建設工程學部,遼寧 大連 116024；2.大連理工大學 海岸和近海工程國家重點實驗室，遼寧 大連 116024;3.深圳中廣核工程設計有限公司，廣東 深圳 518031; 4.核電安全監控技術與裝備國家重點實驗室，廣東 深圳 518172)

能源是支撐當今社會經濟發展的命脈，當今中國經濟迅速發展,能源需求旺盛。傳統化石能源過量開采和使用導致生態環境日益惡化。核電作為一種清潔能源在優化當今能源結構，改善大氣環境方面起到重要作用。防止核泄露是核電廠設計中首要關注的問題。安全殼作為重要核安全相關結構是防止放射性物質泄露的最后一道屏障[1]，設計安全殼結構時需要考慮諸如地震、龍卷風等外部事件以及內部失水事故等工況[2]。安全殼在設計基準工況和超設計基準工況下的完整性可以采用極限內壓承載力進行量化[3-5]。目前對于安全殼在事故工況下內壓承載力研究要集中在試驗研究和數值分析。

關于安全殼試驗研究方面，Hessheime等[6]在桑迪亞實驗室展開了1∶4預應力混凝土安全殼模型試驗。Parmar等[7]進行了1∶4 BARCOM安全殼模型試驗，并且研究了安全殼在不同內壓下的力學響應。Twidale等[8]開展了1∶10安全殼模型試驗，并校核安全殼設計規范。Rizkalla等[9]測試了加拿大某核電廠1∶14安全殼模型。趙樹明等[10]對秦山核電二期安全殼進行了整體性試驗，研究了在安全殼在0.35 MPa內壓和145 ℃極端工況下的力學性能。

由于安全殼模型試驗成本高并且周期長，采用數值模擬方法模擬安全殼非線性行為通常是不可避免的[11]。Hora等[12]提出了安全殼結構長期行為的數值模擬方法。Hu等[13]在考慮長期預應力損失的情況下采用有限元分析了安全殼的極限承載力和破壞模式。 Lee等[14]在考慮混凝土受拉剛化效應的基礎上評估了安全殼極限內壓承載力。Saudy[15]采用有限元方法分析了CANDU堆型安全殼的極限承載力。Lundqvist等[16]評價了安全殼各種預應力損失計算方法，并提出改進預測方法。Shokoohfar等[17]采用塑性損傷模型進行安全殼在高溫和內壓共同作用下的力學性能研究。Tong等[18]采用ANSYS軟件詳細模擬了在預應力和溫度作用下安全殼內壓承載力。Huang等[19]采用有限元軟件真實模擬了安全殼在使用荷載以及假定事故工況下安全殼的性能并且考慮了諸如徐變，收縮等時變因素的影響。Zhang等[20]采用有限元軟件模擬5種蒸汽爆發工況下安全殼結構的完整性。

上述研究均屬于安全殼性能確定性研究，這些研究可以確定安全殼的極限承載力和相應的破壞模式，但是這些確定性試驗和數值模擬研究無法識別安全殼潛在的風險，以及無法了解安全殼在事故工況下的可靠性。由于安全殼結構在設計階段、材料制備及施工和交付使用過程中不可避免的存在各種不確定性因素，從而導致安全殼結構的承載力具有很大的不確定性。同時傳統安全殼設計規范ASME2007[21]采用容許應力法進行設計，通過限制安全殼中鋼筋、混凝土以及預應力筋的應力來保證安全殼的具有一定安全裕度。這種設計方法無法量化安全殼在各種超設計基準工況下的可靠性。同時目前國內外采用概率風險方法對安全殼在事故工況下性能評估的相關報道較少。Ellingwood等[22-25]提出將可靠度理論用于核電廠結構分析。Pandey[26]根據安全殼在役監測數據采用可靠度的方法分析了有粘結預應力安全殼的完整性。Kim等[27]評價了長期預應力損失下安全殼的可靠性，并且分析了安全殼不同位置處的可靠性。上述關于安全殼可靠性研究是以安全殼結構功能函數為簡單顯示函數為前提。由于在實際情況下安全殼的功能函數十分復雜且無法顯示表達。此外關于預應力混凝土安全殼在事故工況下性能評估的相關報道很少。

為此本文采用拉丁超立方抽樣技術對安全殼在事故工況下的易損性進行評估，考慮了材料不確定性下安全殼發生功能性失效、結構性失效模式下的易損性，同時分析不同安全殼變截面位置處發生功能性失效時的易損性。最后采用概率評估方法分析了安全殼在嚴重事故工況下的性能。

1 安全殼有限元模型

1.1 安全殼幾何模型及預應力體系

本文分析的某核電廠預應力混凝土安全殼由穹頂、環梁、扶壁、筒體以及底板等主要部分組成。同時為了保證安全殼的密封性在安全殼內表面設置6 mm鋼襯里。安全殼筒體部分內徑為22.5 m，壁厚為1.2 m。穹頂由曲率半徑為30 m和8 m這2段圓弧構成，穹頂壁厚為1.0 m。安全殼筒體上有設備閘門、普通人員閘門和應急人員閘門3個主要洞口。在筒體的豎向和環向布置了2層普通鋼筋和預應力筋。

安全殼的預應力體系由水平預應力鋼束、純豎向鋼束以及由部分豎向預應力沿著穹頂延伸的倒U型鋼束組成。安全殼的預應力體系采用后張法施工。水平預應力鋼束采用兩端張拉，并錨固在環梁上。豎向預應力鋼束采用一端張拉，一端錨固在環梁上，另一端錨固在基礎底板的預應力筋廊道上。倒U型鋼束也采用一端張拉，一端錨固在對邊的環梁上，另一端也錨固在預應力筋廊道上。安全殼結構簡圖如圖1所示。

圖1 安全殼幾何模型示意Fig.1 Simplified diagram of containment geometry

1.2 預應力損失計算及模擬方法

核電廠安全殼的預應力損失會對安全殼產生不利影響，因此有限元分析需要考慮預應力損失。核電廠預應力混凝土安全殼預應力體系采用后張法施工，預應力損失由長期損失和短期損失2個部分構成。短期預應力損失包含混凝土彈性壓縮、摩擦損失、預應力筋錨具變形造成的損失；長期損失包含徐變、收縮、松弛造成的損失。本文采用實體力筋法模擬混凝土預應力。預應力筋向混凝土的傳遞通過降溫法模擬，預應力筋的有效預應力施加在預應力筋的每一個節點上。預應力筋每個節點的降溫值為：

ΔT=σpe/(αEp) (1)

式中：σpe為考慮短期損失和長期損失后預應力筋的有效預應力；α為預應力筋的線膨脹系數，一般取為10-5/℃；Ep為預應力筋彈性模量。

1.3 材料本構關系模型

1.3.1 混凝土本構關系

(3)

式中：σt、σc分別表示混凝土的拉應力和壓應力；dt和dc分別表示混凝土受拉和受壓損傷因子，具體計算方法見文獻[28]；E0表示混凝的初始彈性模量；εt和εc分別表述混凝土拉應變和壓應變。混凝土單軸受壓本構關系參考文獻[29]。混凝土受拉本構關系采用文獻[30]建議的模型。混凝土單軸受壓本構關系如圖2所示。

圖2 混凝土單軸受壓本構關系Fig.2 Uniaxial compression stress-strain curve of concrete

1.3.2 鋼材本構關系

普通鋼筋、鋼襯里、電氣貫穿件均采用理想彈塑性本構關系模型[31-32]。預應力筋本構關系采用雙線性本構關系。同時普通鋼筋、鋼襯里、電氣貫穿件和預應力筋的泊松比均取為0.3。非預應力鋼材本構關系和預應力筋本構關系如圖3所示(其中fy、εy為非預應力鋼材屈服強度和屈服應變；Ey為非預應力鋼材彈性模量；fp、fu為預應力筋的屈服強度和極限強度；Ep為預應力筋的硬化模量)。

圖3 鋼材本構關系Fig.3 Constitution relationship of steel

1.4 有限元網格劃分、接觸及網格敏感性分析

安全殼混凝土采用以C3D8R單元為主，洞口少量采用C3D6單元劃分網格。預應力筋采用T3D2單元，鋼襯里和電氣貫穿件采用以S4R單元為主，洞口位置采用S3類型單元劃分網格。普通鋼筋采用SFM3D4和SFM3D3單元來劃分網格。所用鋼筋層間距、方位均根據安全殼實際施工圖確定。普通鋼筋單元和預應力筋單元全部嵌入到混凝土單元中，不考慮它們與混凝土的相對滑移。鋼襯里采用ABAQUS軟件中蒙皮的方法來模擬，鋼襯里與混凝土共用節點且兩者沒有相對滑移[33]。安全殼各個部件網格劃分如圖4所示。由于混凝土開裂特性，混凝土結構有限元計算結果對網格尺寸敏感，為此本文進行了安全殼網格尺寸敏感性分析。圖5給出了3種網格尺寸下安全殼最大位移響應隨內壓變化的曲線。從圖5來看，網格尺寸0.6 m與網格尺寸0.8 m計算結果相差很小，網格尺寸1.0 m計算結果與網格尺寸0.6 m和0.8 m在內壓小于0.8 MPa時，計算結果相差甚小，但在內壓超過0.8 MPa以后差異明顯。因此，本文選擇0.8 m網格尺寸進行有限元分析。

圖5 不同網格尺寸下安全殼最大位移與內壓關系曲線Fig.5 Maximum displacement versus internal pressure curves under different mesh size

1.5 隨機變量統計特性

本文選取最為常見的材料不確定性作為不確定性研究對象，研究在考慮材料不確定性下安全殼的內壓易損性。 由于混凝土彈性模量、抗拉強度與抗壓強度相互關聯，因此本文在進行易損性分析時考慮了它們的關聯關系。混凝土彈性模量Ec、抗拉強度ft與抗壓強度fc之間關系可以參考歐洲規范計算。為了簡化分析，在本文有限元分析中材料密度和泊松比不作為隨機變量。所有隨機變量概率統計特性列在表1中。

表1 隨機變量概率統計特性Table 1 Statistical characteristic of random variables

1.6 安全殼破壞準則

安全殼的破壞模式可以分為功能性失效和結構性失效2種，安全殼的功能性失效可以定義為安全殼失去防止泄露的功能。結構性失效是從結構極限強度方面定義安全殼的破壞狀態，結構性失效可以認為是結構的整體變形、應變超過相應的限值[34]。本文采用文獻[35]建議的預應力混凝土安全殼破壞準則。具體破壞準則可以歸納為：

1) 鋼襯里應變達到0.4%(安全殼功能性失效準則)；

2) 普通鋼筋應變達到0.4%(安全殼結構性失效準則)；

3) 預應力筋總應變達到0.8%(安全殼結構性失效準則)。

2 有限元計算結果分析

2.1 自重和預應力工況下安全殼位移響應

圖6給出了安全殼在施加重力荷載和預應力工況下安全殼整體位移分布云圖。在施加重力荷載和預應力工況后，安全殼筒體部位以向內收縮變形為主，穹頂以向外膨脹變形為主。此時安全殼穹頂的位移最大，穹頂產生的向外變形約為18.6 mm，安全殼向內收縮變形大小約為10.8 mm。

圖6 安全殼預應力加重力荷載工況下位移分布Fig.6 Displacement distribution of containment for prestressing and gravity condition

2.2 安全殼內壓作用下非線性位移響應

安全殼對應混凝土剛發生開裂時，最大位移發生在設備閘門兩側小范圍區域，此時最大位移約為20.24 mm，安全殼其他部位的位移相對較小(如圖7(a)所示)。鋼襯里達到破壞準則時，安全殼位移發生重分布，同時安全殼整體位移出現大幅度增長，此時安全殼的最大位移發生在穹頂頂部，安全殼最大位移約為120.8 mm(如圖7(b)所示)。在普通鋼筋和預應力筋達到相應破壞準則時，安全殼位移分布規律與鋼襯里達到相應破壞準則時相同。但在普通鋼筋達到破壞準則時，安全殼的最大位移較鋼襯里達到相應破壞準則時有較大幅度的增長(如圖7(c)所示)。預應力筋達到相應破壞準則時，安全殼整體最大位移與普通鋼筋達到相應破壞準則時對應的最大位移增長幅度較小(如圖7(d)所示)。

圖7 安全殼不同受力狀態下位移分布Fig.7 Displacement distribution of containment under different loading conditions

3 安全殼易損性分析

3.1 拉丁超立方抽樣技術

Monte-Carlo作為一種經典的不確定性分析方法，廣泛運用在工程領域。由于傳統Monte-Carlo方法需要計算工作量大，計算效率低。尤其對于安全殼這種復雜結構采用傳統Monte-Carlo方法進行分析計算工作量巨大。為此本文采用高效率Monte-Carlo抽樣方法進行易損性分析。拉丁超立方抽樣技術是一種分層抽樣方法，可以有效避免重復抽樣，能以較小的樣本量反映總體的變異性。因此，采用拉丁超立方抽樣技術可以大大減少樣本計算數量，節約計算時間。分層抽樣法與重要抽樣方法相似，它們都是使貢獻大的抽樣更多的出現。但是分層抽樣方法不改變原來的概率密度函數，而是將抽樣區間分成一些子區間，并使各個子區間的抽樣點數不同，對貢獻大的子區間抽取更多的樣本[36]。拉丁超立方抽樣方法主要分為以下2個步驟：

(4)

圖8 拉丁超立方抽樣示意Fig.8 Illustration of Latin Hypercube Sampling method

3.2 安全殼易損性計算結果

根據上述破壞準則可以得到安全殼相應的內壓承載力。圖9給出了安全殼發生功能性失效、結構性失效時以及不同變截面位置發生功能性失效對應的內壓承載力分布。

圖9 安全殼內壓承載力分布Fig.9 Pressure capacity distribution of containment

文獻[37]認為安全殼的內壓承載力服從對數正態分布，因此安全殼在事故工況下的內壓易損性可以表達為：

(5)

式中：p為安全殼的內壓荷載；pm為安全殼內壓承載力均值；β為對數標準差。根據上述安全殼內壓承載力服從對數正態分布假定，本文采用最大似然法得到安全殼功能性失效、結構性失效以及不同變截面位置功能性失效準則下對應易損性曲線統計參數。表2、表3分別給出了安全殼發生功能性失效、結構性失效以及安全殼不同變截面位置發生功能性時對應易損性曲線統計參數。

表2 安全殼功能性、結構性失效易損性統計參數Table 2 Statistical parameters of fragility for containment functional failure and structural failure

表3 安全殼不同變截面位置功能性失效易損性統計參數Table 3 Statistical parameters of fragility for different section transition zone functional failure

同時，不同置信水平的安全殼內壓承載力計算公式為：

pQ=pm·eΦ-1(Q)β

(6)

式中：pQ為不同分位值的內壓承載力；Q為分位值水平；Φ-1(·)表示正態分布的反函數。

圖10、表4和表5分別給出安全殼功能性失效、結構性失效以及安全殼不同變截面位置處功能性失效對應的易損性曲線以及相應5%分位、95%分位對應的內壓承載力。從圖10(a)和表4綜合來看，安全殼的破壞失效由功能性失效控制。安全殼發生功能失效對應的內壓承載力與發生結構失效對應的內壓承載力差異較大。基于普通鋼筋失效準則與預應力筋失效準則對應的內壓承載力相差較小。從圖10(b)和表5結合來看，不同變截面位置處發生功能性失效的次序為環梁底部、環梁頂部、穹頂變截面位置、截椎體部位。變截面部位發生功能性失效比安全殼整體發生功能性失效時的內壓承載力大，安全殼環梁底部、穹頂變截面位置、截錐體部位發生功能失效對應的內壓承載力大于安全殼整體發生功能性失效對應的內壓承載力。(同時本文偏于安全，取安全殼整體發生功能性失效時對應5%分位和95%分位的內壓承載力分別作為安全殼的上限和下限內壓承載力)。因此可以得到安全殼的下限和上限內壓承載力分別為1.234 9 MPa和1.362 6 MPa。

圖10 安全殼易損性曲線Fig.10 Fragility curves of containment

表4 安全殼功能性失效、結構性失效對應5%、95%分位內壓承載力Table 4 Pressure capacity of fragility for containment functional failure and structural failure with 5%,95% confidence level MPa

表5 安全殼不同變截面位置功能性失效5%、95%分位內壓承載力Table 5 Pressure capacity of fragility for different section transition zones with 5%,95% percenticle MPa

3.3 安全殼嚴重事故工況下性能評估

目前關于安全殼嚴重事故工況下的性能指標評價的相關研究報道較少。文獻[39]指出安全殼確定性的性能指標為在堆芯熔化24 h后安全殼仍然能保持其完整性和安全殼在嚴重事故工況下的條件失效概率(CCFP)不超過10%。安全殼在事故工況下的條件失效概率可以表達為：

(9)

式中：Ffail(p)為安全殼的內壓易損性曲線；f(p)為安全殼事故壓力分布曲線。

由于本文研究安全殼的設計壓力與文獻[40]中設計壓力不同，因此本文采用修正文獻[40]中給出的事故內壓分布曲線分布的參數得到適合本文安全殼的事故內壓分布易損性曲線。修正后適合本文使用的安全殼事故壓力服從對數正態分布，事故壓力均值為0.663 7 MPa，對數標準差為0.3。本文偏于安全的采用3.2節安全殼整體發生功能性失效對應的易損性曲線作為安全殼嚴重事故工況下的內壓易損性曲線，對卷積積分格式的式(9)進行數值計算得到安全殼的條件失效概率為0.01，說明本文分析的安全殼可以很好地滿足事故工況下的性能指標。

4 結論

1)安全殼破壞由功能性失效控制，且對應安全殼發生功能性失效時對應的內壓承載力與發生結構性失效時對應的內壓承載力相差較大，就結構性失效而言，普通鋼筋失效準則、預應力筋失效準則對應的內壓承載力相差較小。

2)安全殼不同變截面位置發生功能失效的次序為環梁底部、環梁頂部、穹頂變截面位置、截椎體部位，環梁底部、環梁頂部、穹頂變截面位置、截椎體部位發生功能性失效時對應的內壓承載力比安全殼整體發生功能性失效對應的內壓承載力高，說明這些部位相對安全殼整體性失效是偏于安全的。

3)安全殼在嚴重事故工況下的條件失效概率為0.01遠小于性能指標0.1，說明該安全殼可以很好的滿足嚴重事故工況下要求的性能要求。