大型商用飛機撞擊核電站屏蔽廠房荷載研究

劉晶波， 鄭文凱

(清華大學 土木工程系，北京 100084)

至2011年底我國核電機組已建15臺、在建26臺，擬建21臺[1]。在核電快速發展的同時，核電的安全性亦愈受重視。鑒于核電站的特殊性及核燃料泄漏后果的嚴重性，須重點考慮核電站安全性，尤其須確保核電站受火災、地震、海嘯、飛機撞擊等極端情況的安全[2]。

“9·11”巨大災難引起公眾對大型商用飛機撞擊重要構筑物的關注。隨核電建設發展及大型商用飛機數量劇增，飛機撞擊核電站問題亦引起重視。核電站設計中已開始考慮飛機撞擊屏蔽廠房因素[3]。對此可直接用飛機有限元模型撞擊核電站屏蔽廠房有限元模型方法進行分析計算，該方法除需建立核電站結構模型外，還需建立復雜的非線性飛機模型、考慮飛機與核電站動接觸關系及兩者的動力相互作用，計算較繁復。工程設計中希望能得到飛機撞擊荷載時程曲線，利用該時程曲線直接加載考慮飛機的撞擊作用，方便實用。但目前我國對飛機撞擊荷載研究較少[4-5]。

本文利用ANSYS/LS-DYNA建立非線性飛機與核電站屏蔽廠房的有限元模型模擬撞擊過程，考慮非線性、大變形、動接觸等因素，通過發動機撞擊鋼筋混凝土板試驗驗證材料的本構模型及參數，對飛機網格分析確定其尺寸，對飛機以200 m/s速度撞擊核電站進行分析，獲得飛機撞擊作用下核電站變形、破壞形式等規律。對假設的核電站采用線彈性、剛性本構模型，分析、對比核電站采用不同本構模型時撞擊荷載時程，結果顯示核電站剛度相對飛機較大，可采用撞擊剛性體方法[6]確定飛機撞擊合力的荷載時程曲線。但該方法無法給出荷載分布形式，僅建議將荷載均勻分布在機身作用面上。分析表明該分布方式結果與實際撞擊結果相差較大。結論對確定核電站設計中飛機撞擊荷載時程曲線與分布具有實用參考價值。

1 材料本構模型

大型商用飛機與核電站屏蔽廠房發生撞擊時作用時間短、峰值大、應變率高，各種材料的力學性能與靜載作用相比差距較大，故確定各種材料的本構模型較關鍵。大型商用飛機主要由鋁合金與鋼材組成。核電站屏蔽廠房主要材料為鋼材及混凝土。本文對鋁合金、鋼材、混凝土材料的本構模型及參數取值進行分析研究。

1.1 鋼材及鋁合金材料本構

或飛機或核電站結構中的金屬材料碰撞時材料的應變率效應十分明顯，故選用考慮屈服、強化及應變率效應的Cowper-Symonds本構模型[7]，該模型為碰撞分析中常用的金屬本構，屈服函數為：

(1)

1.2 混凝土材料本構

由于核電站主要受力材料為混凝土，撞擊荷載下需研究混凝土非線性本構關系及損傷，并考慮應變率效應。HJC[8]、Brittle_Damage[9]及Concrete_Damage[10]為動力分析的混凝土本構。本文選考慮鋼筋作用、應變率效應、損傷效應、應變強化及軟化作用的Concrete Damage本構模擬混凝土。該本構專為分析混凝土受爆炸及撞擊作用，并定義三獨立強度面：

(2)

其中：i分別為y,m,r代表初始屈服面、極限強度面、殘余強度面；p=-I1/3為靜水壓力，可通過狀態方程定義p與體積應變之關系；aji(j=0,1,2)為常數，由實驗確定。

后繼屈服面可利用參數η通過線性內插獲得，參數η可反映屈服面強化、軟化，且為與累計有效塑性應變相關函數，在0～1間變化。達到初始屈服面、未達極限強度面通過式(3)內插計算后繼屈服面，此時η由0變到1；達到極限強度面后開始軟化，采用式(4)內插后繼屈服面，此時η由1變到0。

Δσ=ηΔσm+(1-η)Δσy

(3)

Δσ=ηΔσm+(1-η)Δσr

(4)

本模型參數較多，需較多實驗確定，本文取實驗參考值[7,11]。

1.3 材料模型及參數驗證

1.3.1 驗證試驗

文獻[12]用簡化發動機分別以100 m/s、150 m/s、215 m/s速度對正方形鋼筋混凝土板撞擊試驗。正方形板邊長1.5 m，厚12 cm，板中縱筋配筋率0.4%，間距60 mm，直徑6 mm，無箍筋。發動機見圖1，總重3.5 kg。本文擬利用此試驗對模型中混凝土、鋼材本構及相關參數進行驗證。

圖1 試驗用發動機模型(單位：mm)

試驗測得混凝土圓柱體抗壓強度23.5 MPa，密度2300 kg/m3；板中縱筋屈服強度447.2 MPa，彈性模量2.05×105MPa；發動機鋼材屈服強度411.9 MPa，彈性模量2.14×105MPa。有限元模型中混凝土材料采用CONCRETE_DAMAGE本構，本構中應變率效應曲線用LS-DYNA經大量實驗確定的推薦曲線[7]。鋼筋、鋼板均采用Cowper-Symonds本構，參數見表1。為模擬混凝土材料失效，用單元失效的EROSION算法，該算法可通過定義應力、應變等多種閥值控制單元失效，計算中需根據計算結果及試驗結果對照確定相應閥值，由于混凝土應力變化大且本構曲線存在下降段，故選最大主應變控制單元失效，通過試算取最大失效主應變為0.17。

表1 鋼板及鋼筋材料參數

據試驗建立有限元模型見圖2?；炷劣脝卧叽?0 mm的Solid164實體單元離散；鋼筋用Link160桿單元離散，桿單元尺寸30 mm。鋼筋、混凝土用共用節點建模。發動機鋼板用Shell163殼單元離散，長度方向單元尺寸10 mm。模型單元總數285872個，鋼筋混凝土板沖擊背面四角點處，在200×200 mm范圍內用固支約束。

用罰函數方法處理接觸界面以模擬發動機對鋼筋混凝土板的撞擊作用，該方法根據接觸時主、從界面相關節點的穿透情況引入與穿透深度、剛度成正比的界面接觸力以保證接觸界面不發生穿透。

在建立發動機與鋼筋混凝土板有限元模型后，分別定義發動機與鋼筋混凝土板中鋼筋、混凝土的接觸方式。發動機與鋼筋接觸方式用發動機為主界面的自動點面接觸算法；發動機與混凝土接觸方式用自動面面接觸算法。罰函數方法中罰函數因子據實際接觸中穿透情況調節，摩擦系數據實際接觸表面情況確定。本文取罰函數因子1.0以防止穿透，摩擦等其它系數用默認值。

1.3.2 試驗與數值模擬分析

圖3(a)為數值模擬發動機以100 m/s速度撞擊時混凝土板正面被撞擊成坑，背面未破壞，與圖3(b)試驗照片形態大小一致。圖4為發動機以150 m/s速度撞擊時板正面、背面的破壞情況，正面混凝土破壞成坑，背面混凝土破壞脫落縱筋外露。圖5為215 m/s撞擊時穿透混凝土板，撞擊處鋼筋斷裂，混凝土破壞，數值模擬發動機穿透后的剩余速度約52 m/s，與實驗所測55 m/s較接近，破壞孔洞形態與實驗一致。

圖2 發動機撞擊鋼筋混凝土板有限元模型

圖4 150 m/s撞擊時混凝土正面、背面破壞

撞擊試驗及數值模擬結果對比見表2。由表2看出，鋼筋混凝土板的破壞類型、形態與實驗一致，表明所用材料模型能較好模擬鋼材、混凝土被沖擊作用性質，且材料失效判斷類型及參數取值合理。

表2 試驗與數值模擬結果對比

2 核電站與飛機模型

2.1 核電站模型

以某型號核電站為原型并結合國際常用核電站結構尺寸建立核電站屏蔽廠房模型。筒體結構外半徑22 m，高44 m，頂部近似半球形外半徑22 m，屏蔽廠房由1 m厚混凝土雙面包13 mm厚鋼板組成。結構總高66 m，混凝土用Solid164實體單元劃分，鋼板用Shell163殼單元劃分，考慮計算成本，單元邊長25 cm，單元總數1018368，其中實體單元數678912。結合核電站屏蔽廠房工程實際，鋼板與混凝土間認為連接可靠，用共用節點建模忽略兩者間滑移，結構底部為固定邊界?；炷敛牧媳緲嫴捎肅ONCRETE_DAMAGE本構，混凝土單軸抗壓強度48 MPa，鋼板用Cowper-Symonds本構，其它參數同前。

2.2 飛機模型

選飛機模型關鍵控制量為總質量、發動機質量及質量分布。常選撞擊飛機重量、飛行量應能包絡大多數的商用飛機[13]。由于難以獲取精確的飛機構造等數據，本文以Boeing 767-200ER大型商用飛機尺寸、重量為參考[14]。飛機重80 t，發動機重約10 t，飛機總長48.5 m，添加地板及橫隔板模擬內部結構。飛機全部采用殼單元離散；針對EROSION算法在單元失效時刪除單元，飛機網格尺寸對撞擊結果影響較大，故進行網格尺寸分析。圖6為網格尺寸分別為1 m，0.25 m，0.125 m時飛機撞擊剛體的撞擊力合力，隨網格尺寸減小計算趨于穩定。本文選飛機網格尺寸0.25 m，計算結果穩定且計算成本較小，飛機單元總數9148，飛機外形尺寸及有限元模型見圖7。

圖6 網格尺寸對撞擊力影響分析

3 飛機撞擊核電站計算模擬

飛機撞擊核電站屏蔽廠房分析中涉及到飛機的撞擊位置、角度及速度。文獻[15]指認為最不利撞擊位置在筒體與半球交接處,而垂直撞擊為屏蔽廠房整體分析中最不利作用方式，撞擊速度取200 m/s。本文計算中飛機以200 m/s速度垂直撞擊核電站筒體與半球交接處，見圖8。

用上述有限元模型，撞擊模型單元總數1027516，撞擊中采用部件組定義接觸為自動面面接觸，其中罰函數作用系數為1.0以保證合理接觸，材料模型及參數同上。本文僅考慮飛機撞擊作用，未考慮飛機燃油燃燒及爆炸作用。撞擊過程見圖9，撞擊結果顯示飛機機身撞擊部分已完全破壞，核電站模型撞擊部位變形明顯，機身及發動機巨大撞擊導致被撞擊區域混凝土破壞，發動機在飛機撞擊中影響較大，其它區域反應較

小。圖10為核電站結構距被撞擊中心點上方3 m、5 m、7 m、9 m四個點及頂點的位移曲線，隨距離撞擊中心點越遠撞擊影響迅速減小，核電站頂點位移最大值僅15 mm。飛機機身未破壞部分速度見圖11，飛機撞擊力合力見圖12，最大撞擊力合力約170 MN。

圖9 飛機撞擊過程

圖10 結構位移曲線

4 撞擊荷載分析與討論

影響飛機撞擊力的因素較多，如飛機質量大小、質量分布及材料性質、核電站屏蔽廠房外形尺寸及材料性質、飛機速度等。為進一步分析撞擊荷載作用，用非線性飛機模型分別撞擊剛性、線彈性核電站模型，獲得撞擊力曲線見圖13。由圖13看出，三種情況整體作用時間、形狀及峰值相近，說明與飛機相比，核電站剛度較大，可設為剛性體或彈性體估算撞擊合力。

由于核電站剛度相對飛機大，可采用撞擊剛性體方法[6]分析飛機撞擊合力的荷載時程曲線。該方法主要假設：飛機垂直撞擊剛性體目標，飛機為質量沿長度方向分布的1維結構。據動量定理導出飛機撞擊力計算式為：

F(t)=Pc(x)+μ(x)(dx/dx)2

(5)

其中：F(t)為飛機撞擊力；x(t)為機頭至t時刻壓屈作用點距離；Pc(x)為飛機在x處壓碎力；μ(x)為在x處單位長度飛機質量。用文獻[6]方法求解飛機初速度200 m/s時的撞擊荷載時程見圖14。由圖14看出，利用Riera方法所求撞擊合力的荷載時程曲線與彈塑性數值模擬所得撞擊合力除第2峰值列處相差較大外總體形狀吻合較好，故可用該方法計算撞擊合力的荷載時程曲線。但Riera方法無法給出撞擊荷載分布形式，僅將荷載平均分布到飛機機身作用面積內，按均勻分布給核電站加載。模擬結果顯示采用此加載方式計算時均布荷載邊緣混凝土先遭破壞，破壞形態及原因與飛機撞擊非線性核電站模擬結果差別較大。圖15為均布加載時核電站結構上荷載作用中心點上方3 m、5 m、7 m、9 m四點與頂點位移曲線，與圖10相比看出兩者差異較大，撞擊處位移及整體反應差別較大，說明不能簡單將荷載均勻分布在機身作用面積上，需考慮不同部位各自荷載大小及分布。

圖13 剛性體及線彈性體撞擊力

5 結 論

(1)本文通過發動機撞擊鋼筋混凝土板試驗確定CONCRETE_DAMAGE與Cowper-Symonds本構能較好模擬撞擊作用下混凝土與金屬材料性質，并確定接觸算法及相關參數。通過用不同網格尺寸飛機撞擊剛性體分析確定飛機模型合理網格尺寸。

(2)通過對大型商用飛機以200 m/s速度撞擊鋼板混凝土核電站屏蔽廠房的非線性動力時程分析結果顯示，核電站被撞擊部位變形明顯，機身及發動機撞擊發生破壞，其它區域變形較小，最大撞擊合力約170 MN。

(3)將核電站用彈塑性模型時撞擊合力荷載時程曲線與用剛性、線彈性體核電站模型撞擊荷載對比，三種情況下曲線相近說明，與飛機相比，核電站剛度較大，可設核電站為剛性體或彈性體計算飛機撞擊力合力，可用Riera方法推導出總撞擊力荷載時程曲線。而該方法無法確定荷載分布形式，所用在機身作用面均勻加載方式所得結果與飛機撞擊核電站非線性模型結果差別較大，故需對荷載分布形式進一步研究。

參 考 文 獻

[1]趙小輝,鄒樹梁,劉 永.內陸核電發展形勢分析[J]. 南華大學學報,2012,13(3):1-7.

ZHAO Xiao-hui,ZOU Shu-liang,LIU Yong. Analysis of inland nuclear power development situation[J]. Journal of University of South China, 2012, 13(3):1-7.

[2]王玉嵐.核電站安全殼在外來事件作用下的安全防護性能研究進展[J]. 電網與清潔能源,2009,25(10):74-77.

WANG Yu-lan. Research progress of security performance of nuclear power plant suffering external events[J]. Power System and Clean Energy, 2009, 25(10):74-77.

[3]湯 搏.關于核電廠防大型商用飛機撞擊的要求[J].核安全, 2010,3:1-16.

TANG Bo.Requirements for nuclear power plants against large commercial aircraft impact[J].Nuclear safety,2010,3:1-16.

[4]左家紅.秦山核電廠安全殼在飛機撞擊下的非線性分析[J].核科學與工程,1992,12(1):35-42.

ZUO Jia-hong. Nonlinear analysis of Qinshan nuclear power plant containment under impacting of the aircraft[J]. Nuclear Science and Engineering, 1992,12(1):35-42.

[5]王遠功,余愛萍. 飛機撞擊核反應堆安全殼荷載-時間曲線的確定[J]. 核科學與工程,1991,11(3):208-215.

WANG Yuan-gong,Yu Ai-ping. The determination of the load of the aircraft hitting the reactor containment[J]. Nuclear Science and Engineering, 1991, 11(3):208-215.

[6]Riera J D. On the stress analysis of structures subjected to aircraft impact forces[J].Nuclear Engineering and Design, 1968, 8(4):415-426.

[7]LS-DYNA Keyword User’s Manual,Version971,2007．

[8]Holmquist T J, Johnson G R. A computational constitutive model for concrete subjected to large strains, high strain rates, and high pressures[R]. Canada, Proceedings of the fourteenthInternational Symposium on ballistic, 1993.

[9]Govindjee S, Kay G J,Simo J C. Anisotropic modelling and numerical simulation of brittle damage in concrete[J]. International Journal for Numerical Methods in Engineering,1995,38(21):3611-3633.

[10]Malvar L J, Crawford J E, Wesevich J W,et al.A plasticity concrete materialmodel for DYNA3D[J].International Journal of Impact Engineering, 1997, 19( 9-10):847-873.

[11]Crawford J E, Magallanes J M, Lan S,et al. User’s manual and documentation for release III of the K&C concrete material model in LS-DYNA, TR-11-36-1[R].Technical Report, Karagozian & Case, Burbank, CA, 2011.

[12]Muto K, Tackihawa H, Sugano T, et al. Experimental studies on local damage of reinforced concretestructures by the impact of deformable missiles, part 1: outline of test program and small-scale tests SmiRT 10[R].Los Angeles, California: American Association For Structural Mechanics in Nuclear Reactors, 1989.

[13]EPRI.Deterring terrorism: aircraft crash impact analyses demonstrate nuclear power plant's structural strength [R]. America: Electric Power Research Institute, 2002.

[14]Boeing Commercial Airplanes[C]. http://www.boeing.com

[15]Prabhakar G, Ranjan R, Paul M K,et al.Analyses of aircraft impact on containment structure[C]. 5thAsia-Pacific Conference on Shock & Impact Loads on Structures. Changsha China, 2003:315-322.