核電廠房安全殼流固耦合動力響應分析

曾 濱，王宇威，潘鉆峰，許 慶，薛 衛(wèi)

(1.中冶建筑研究總院有限公司，北京 100088；2.同濟大學土木工程學院，上海 200092；3.華龍國際核電技術有限公司，北京 100036)

非能動安全殼水儲存箱內水體在地震荷載激勵下會產生晃動，水體與水箱壁的交互作用被稱之為流固耦合，這種效應會對安全殼整體結構的動力響應產生影響。

目前，不少學者對大型儲液結構的流固耦合問題開展了研究，簡化附加質量法最早是由Housner[1]提出，Housner 模型常用在水槽結構和各類儲液罐結構的流固耦合效應分析中[2-3]，AP1000 內置換料水箱即采用了Housner 模型進行地震作用動力分析[4]，這種方法在實際的運用過程中，由于計算的復雜，往往會進一步簡化忽略水體的晃動對流作用，僅考慮脈沖質量，簡化附加質量法是一種偏向于保守的計算方法。Donea 等[5]、Souli 和Zolesio[6]提出了任意拉格朗日-歐拉耦合(ALE)法，寶鑫等[7]利用此方法結合有限元研究了土-結構動力相互作用對儲液結構地震反應的影響，樓云峰等[8]對流體區(qū)域采用ALE 描述研究了多工況波浪沖擊下浮橋的動力響應，趙春風等[9-12]利用AEL 方法分析了水箱水位高度對核島廠房自振頻率、加速度和位移等動力響應的影響。光滑粒子流體動力學(SPH)法[13-14]是一種新型的、無網格的數值研究方法，SPH 方法的理論基礎是粒子方法，粒子方法是對連續(xù)物理量用多數粒子的集合插值進行求解的數值方法，這些粒子每個都代表了流體的一個微元，并且這些微元都帶有相應的材料屬性(密度、速度、位置等)，然后通過求解整個粒子系統(tǒng)的動力學方程組來模擬流體的流動。SPH 方法作為一種新型的、無網格的數值研究方法，其優(yōu)勢在于粒子系統(tǒng)較為真實模擬水體流動，計算精度高，其不用網格的特點，也不會產生計算溢出的問題，即不會存在因形狀不規(guī)則引起的漏水現象。SPH 方法和有限元的結合[15-16]使得很多實際工程問題[17-19]的研究得到了極大的方便。在流固耦合工程問題研究中也越來越多，Antoci 等[20]用粒子網格對水渠進行模擬并和試驗結果對比，驗證了SPH 法運用于流固耦合模擬的可行性；Groenenboom 和Cartwright[21]用SPH-FE 方法對圓筒跌落于靜水的模擬，觀察到液面的波動情況與真實現象非常相似；Xu 等[22]采用SPH-FE 法分析了AP1000 水箱水位高度對屏蔽廠房動力響應的影響。由于安全殼水儲存箱體型巨大，其水箱的流固耦合動力效應將對結構產生不容忽視的影響，在已有研究基礎上，結合SPH 方法的流固耦合分析有進一步的研究價值。

本文以國內某新型核電廠房安全殼為背景，利用ABAQUS 有限元分析軟件對某新型安全殼進行了建模計算分析，獲取了該型安全殼的前10 階模態(tài)以及水箱不同水位高度工況下的安全殼的前10 階頻率；采用常用自然波Northridge 地震波上海人工合成波作為地震動荷載輸入，對比了簡化附加質量法和SPH 法(分別對應于不考慮、考慮水體晃動效應)安全殼基底剪力、最大相對位移和樓層加速度反應譜等動力響應，并進一步利用SPH 法探究了不同水位高度工況對安全殼的位移、加速度動力響應和樓層反應譜的影響，探尋該型核電廠房安全殼水箱的最合理水位。

1 安全殼整體模型

1.1 主體模型建立

核電廠為單層安全殼預應力混凝土結構，安全殼主要由穹頂、筒身和基礎筏板組成，其中穹頂為半球殼，筒身為圓柱薄殼結構，水箱位于穹頂與筒身相交處，共劃分成3 層，僅上層裝水，其余為設備層，安全殼及水箱幾何特征示意見圖1。安全殼主體標高為-8.240 m～+67.700 m，筒身高51.94 m，穹頂高24.0 m，基礎筏板底部固結，筏板厚4.0 m；筒身內徑45.0 m，外徑48.0 m，厚度1.5 m；水箱標高+37.570 m～+55.500 m，外掛在安全殼外壁，水箱外徑55.9 m，各層高度從上向下依次為5 m、7.5 m 和2.65 m，除儲水層板厚0.8 m，其余均為0.6 m。

圖1 安全殼及水箱示意圖 /mFig.1 Diagrammatic sketch of containment and tank

在建模中，考慮到鋼筋布置的復雜性和避免因接觸帶來的不收斂問題，鋼筋和預應力筋以換算截面積計入了混凝土中。安全殼主體與水箱均采用SR3、SR4 均質殼單元建模，單元尺寸為2000 mm。簡化附加質量法按節(jié)點分擔面積將水體總質量平均成質量單元依附在與水體接觸的水箱壁上，SPH法將水體首先建立成C3D8R 實體單元，然后在指定單元等參方向生成2 個粒子，即每個父單元生成8 個粒子(圖2)，轉換準則基于時間的準則，設置閾值為0。

圖2 SPH 法粒子單元轉換Fig.2 Element conversion of SPH method

1.2 參數分析方案

1.2.1 材料參數

混凝土采用線彈性材料模型，混凝土強度等級C60，彈性模量36 GPa；對于結構的阻尼，對簡化附加質量法和SPH 法均通過定義混凝土的瑞利阻尼系數來設置，瑞利阻尼系數由結構固有頻率逆推獲得，阻尼比按《核電廠抗震設計規(guī)范》(GB 50267-97)[23]取值為0.05。水體采用EOS 狀態(tài)方程，ABAQUS 軟件具有該水動材料模型Us-Up 形式。材料的具體參數見表1，各水位工況混凝土瑞利阻尼系數見表2。

表1 混凝土和水體材料屬性Table 1 Concrete and water material properties

表2 混凝土瑞利阻尼系數Table 2 Rayleigh damping coefficients of concrete

1.2.2 荷載參數

荷載施加設置2 個分析步驟：第1 步對整體模型施加重力荷載g(9.8 N/kg)；第2 步施加地震荷載，地震波選取常用中遠場自然波Northridge 波和遠場上海人工合成波，作三向地震動輸入，并調整峰值加速度至0.3g。圖3 為調整后的Northridge波和上海人工合成波加速度時程曲線。

圖3 地震波加速度時程曲線Fig.3 Acceleration time history

1.2.3 水箱水位高度

儲水箱凈高5 m，設置水位為無水、0.8 m、1.6 m、2.4 m、3.2 m、4.0 m、4.8 m 共7 種工況，工況模型見圖4。

圖4 水位高度工況Fig.4 Water levels

1.2.4 水體單元參數敏感性分析

水體單元的形式，網格尺寸和轉換粒子的數目對計算結果精度有顯著影響，同時網格越密、轉換粒子數目越多，計算效率就越低，水體單元采用粒子形式相比質量塊形式計算時間增長了約20 倍。以Northridge 波作用下工況4 為例，對網格尺寸600 mm、1200 mm 和2400 mm 及每個單元轉換粒子數目8 個和24 個進行對比分析。從基底剪力計算結果來看，網格尺寸600 mm 相比1200 mm 精度提升1.3%，相比2400 mm 精度提升3.7%；轉換粒子數目27 個相比8 個精度提升2.3%，同時計算時間增長了3 倍。兼顧計算效率，網格尺寸選擇1200 mm，粒子轉換數目選擇8 個。

2 模態(tài)和固有頻率分析

提取無水工況下前10 階模態(tài)，由于安全殼為對稱性結構，剔除水平向X、Z方向出現的重復模態(tài)后的前10 階模態(tài)如圖5 所示。第1 階、第2 階表現出水平方向的懸臂梁式彎曲振動；第3 階、第4 階、第6 階、第7 階、第9 階、第10 階為安全殼筒身環(huán)向的多波屈曲振動；第5 階為扭轉振動；第8 階為水箱的豎向伸縮振動。表3 給出了剔除重復模態(tài)后的安全殼前10 階頻率，圖6 是安全殼前10 階固有頻率與水位高度的關系，隨著水位高度的上升，各階模態(tài)振型沒有發(fā)生改變，但各階模態(tài)固有頻率則以不同的速率呈現下降趨勢。其中，第8 階頻率隨水位上升下降的最快，說明水位高度變化還是對水箱本身的影響最為顯著；其次為第1 階、第2 階和第5 階，分別為彎曲和扭轉振動；水位變化對第3 階、第4 階、第6 階、第7 階、第9 階和第10 階的水箱下部筒身的多波屈曲振動則幾乎沒有影響。由此可見，水箱水位的變化對水箱所在位置的動力響應影響最大，結構的整體動力響應次之，而對安全殼下部結構動力響應的影響有限。

表3 各工況的前10 階頻率Table 3 First 10 natural frequencies at various water levels

圖5 安全殼前10 階模態(tài)Fig.5 First 10 modes of containment

圖6 前10 階頻率隨水位高度變化關系Fig.6 First 10 natural frequencies at various water levels

3 簡化附加質量法與SPH 法對比分析

應用SPH 方法充分考慮了水體的晃動效應，同時建立水體簡化質量模型，即不考慮水體的晃動效應，僅將水體作為脈沖質量附著在水箱壁上，通過比較兩種方法下安全殼結構的基底剪力和安全殼上不同高度5 個節(jié)點(圖7)加速度、位移動力響應的差異，研究水箱內水的晃動對安全殼結構的影響。

圖7 節(jié)點位置Fig.7 Positions of monitoring points

3.1 最大基底剪力的比較

利用簡化附加質量法和SPH 法模擬水體計算得到的各水位工況安全殼筏板在兩個地震波作用下最大基底剪力如表4 和表5 所示。從表中的數據可以看出，無論是否考慮水箱水體的晃動效應，最大基底剪力總體上表現出隨著水位高度的增加而增大。各水位工況下，采用SPH 法得到的最大基底剪力均小于用簡化附加質量法所得到的最大基底剪力，這說明SPH 法考慮水的晃動效應后對基底剪力的動力響應是有減弱的，水箱水的晃動對流對安全殼起到了減震耗能的作用。在具體的數值上，Northridge 波計算結果均要大于人工波的計算結果。其中在水位高度為4.0 m 時，Northridge波考慮水體晃動對最大基底剪力的減弱效果最大，達到了9.60%；相同水位下，人工波的計算結果為6.08%。

表4 Northridge 波各水位工況下兩種方法基底剪力比較Table 4 Comparison of base shear under various water levels by two methods

表5 人工波各水位工況下兩種方法基底剪力比較Table 5 Comparison of base shear under various water levels by two methods

3.2 節(jié)點最大相對位移比較

選取水箱上部節(jié)點P1、節(jié)點P2 和節(jié)點P3，對其在Northridge 波和人工波作用下的最大相對位移值比較，最大相對位移指該點位移與基底位移的差值。圖8 是簡化附加質量法和SPH 法模擬的最大相對位移與水位高度的關系圖，從圖8 中可以看出，各位置節(jié)點的最大相對位移曲線形式一樣，簡化附加質量法曲線表現出相對最大位移隨水位高度增大而增大，SPH 法曲線表現為先減小后增大。在各水位工況下，簡化附加質量法模擬結果始終大于SPH 法模擬結果，說明水體的晃動對最大相對位移響應具有降低效果。從圖8 中還可以看出，在高水位時兩種方法的計算差值要比低水位時更大，二者的最大差值出現在水位高度4.0 m 時。

圖8 不同水位工況下兩種方法水平最大位移比較Fig.8 Comparison of maximum displacement under various water levels

3.3 水位高度4.0 m 樓層加速度反應譜比較

綜合最大基底剪力和最大相對位移的比較結果，不難看出SPH 方法的結果要小于簡化附加質量方法的結果，且在水位高度4.0 m 時最為明顯，因此，進一步獲取了水位高度4.0 m 時各節(jié)點在Northridge 波作用下的樓層加速度反應譜，如圖9所示。圖9 給出了P1、P2 和P3 節(jié)點簡化附加質量法和SPH 法在水位高度4.0 m 時的加速度反應譜。從圖9 中可以看出，兩種方法獲取的曲線形式一致，簡化附加質量法的反應譜峰值周期要稍微滯后于SPH 法。利用簡化附加質量法獲取的反應峰值在各個節(jié)點均明顯大于SPH 法的結果，二者的差值隨節(jié)點位置變高而增大。

圖9 水位高度4.0 m 兩種方法各節(jié)點加速度反應譜比較Fig.9 Comparison of acceleration response spectra at water level of 4.0 m

綜合最大基底剪力、最大相對位移和樓層加速度反應譜的簡化附加質量法和SPH 法的比較結果來看，SPH 法考慮水箱內水體的晃動對流作用后，對結構整體的動力響應是具有削弱效果的，即水箱內水體起到了耗能減震的作用。

4 SPH 法水位高度對安全殼動力響應的影響

4.1 水位高度對最大位移響應的影響

圖10 是用SPH 法模擬的各工況下水箱水體自由液面粒子在Northridge 波作用下的水平X向位移云圖，位移云圖可以反映水體的晃動和變形程度。從圖10 中可以看出，工況4～工況6，水位從2.4 m 上升到4.0 m 時，自由液面的水平位移逐漸減??；工況7，水位高度4.8 m 時的位移又略微增大。

安全殼上各節(jié)點的相對基底的最大位移隨水位高度變化的關系如圖11 所示，Northridge 波和人工波作用下的變化規(guī)律大體一致，在數值上Northridge 波更大。從兩幅圖上來看，曲線形式表現為：水位高度0 m～1.6 m，相對最大位移先增大；水位高度1.6 m～4.0 m，相對最大位移呈現出逐漸減小的趨勢；水位高度4.0 m～4.8 m，相對最大位移又出現增大。第一次增大的原因是水體質量增加引起慣性力作用效果要大于水體增加帶來的阻尼效果，中間的下降段則反映了水體晃動的耗能減震作用；第二次增大可能是自由液面高程離水箱頂板過近，影響了水體的自由晃動。從各節(jié)點的位移大小可以看出，安全殼上位置越高，相對最大位移越大，說明安全殼在地震作用下，位置越高，擺動幅度越大。此外，從各節(jié)點相對最大位移的曲線隨水位高度的波動狀態(tài)可以得出：水箱水體晃動對安全殼水箱上部的影響要比安全殼下部的影響要大，具體表現為P1、P2、P3 節(jié)點的曲線波動幅度較大，而底部節(jié)點P4、P5 的曲線形式則相對平緩。

圖11 中安全殼上各節(jié)點的相對最大位移隨水位高度的變化規(guī)律與圖10 水箱水體晃動所展現的現象吻合良好，說明水箱水位在4.0 m 時對相對最大位移動力響應減弱最明顯。

圖10 水箱水體位移云圖Fig.10 Time sequence of displacement for water in water tank

圖11 相對最大位移與水位高度關系Fig.11 Relationship between relative maximum displacement and water levels

4.2 水位高度對最大加速度響應的影響

圖12 是安全殼上各節(jié)點最大加速度和水箱水位高度的變化關系，兩種地震波作用下的最大加速度曲線規(guī)律幾乎一致。整體來看，最大加速度與水位高度的變化關系與最大位移與水位高度的關系略有不同，但在安全殼水箱上部節(jié)點P1、節(jié)點P2 則表現出了的相似的規(guī)律，大致趨勢為先減小后增大，在3.2 m～4.0 m 出現最小值，說明水箱水位在這一區(qū)間時，水體的晃動作用對安全殼上部的加速度動力響應削弱最強。此外，水的晃動效應對安全殼各高度位置最大加速度的影響規(guī)律與相對最大位移相同，表現出水箱上部受到的減震效果要大于水箱下部受到的減震效果。

4.3 樓層反應譜

依據地震波作用下安全殼上各節(jié)點的加速度時程曲線，采用seismosignal 軟件計算了P1、P3 和P5 等布置有相關設備節(jié)點的樓層加速度反應譜，計算結果如圖13 所示。從圖13 中可以看出。當周期約為0.3 s 時，Northridge 波和人工波作用下的安全殼樓層加速度反應譜都出現最大反應，因此設備安裝時應避免設備固有周期在0.3 s 前后一段范圍內。各節(jié)點在各水位工況下的加速度反應譜曲線幾乎一致，除圖13(e)～圖13(f)兩幅圖各水位工況下的反應譜曲線幾乎重疊外，圖12(a)～圖12(d)四圖都能判斷出工況6(水位高度4.0 m)的加速度反應譜最小，而工況1(無水)和工況2(水位高度0.8 m)的反應譜峰值最大。此外，節(jié)點位置越高的反應譜曲線越分化和位置越低反應譜曲線越重疊，也證明了水箱水體的晃動對安全殼上部動力響應的影響要大于下部的結論。

圖12 最大加速度與水位高度關系Fig.12 Relationship between maximum acceleration and water levels

圖13 樓層加速度反應譜Fig.13 Acceleration response spectra

綜上所述，從相對最大位移曲線、最大加速度曲線和樓層加速度反應譜的計算結果可以綜合判斷工況6，水箱水位高度4.0 m 是該新型核電安全殼水箱較為合理的蓄水水位。在這一水位下，水箱水體的晃動作用對安全殼動力響應能起到很好的減震效果。

5 結論

本文研究了國內某新型核電廠房安全殼在Northridge 地震波和上海人工合成波作用下水箱的流固耦合效應對結構動力響應的影響。給出了該型安全殼的前10 階模態(tài)，探究了水箱水位對結構固有頻率的影響；對比了采用簡化附加質量法和SPH 法分別對應不考慮與考慮水箱水的晃動效應對結構動力響應的影響；進一步分析了SPH 法各水位工況與結構動力響應的聯系。主要結論如下：

(1) 水箱水位對結構固有頻率有影響，表明水箱水位變化會影響到安全殼的動力響應，結構前10 階固有頻率隨水位的上升以不同速率降低，其中第8 階水箱豎向伸縮振動的模態(tài)頻率下降速率最快，說明水位影響最直接還是水箱本身，其次為安全殼的整體彎曲和扭轉，對安全殼下部筒身的動力響應則有限。

(2) 對比簡化附加質量法和SPH 法的模擬結果，考慮水體晃動后，基底剪力和最大相對位移響應各水位工況下都出現減弱，基底剪力在4.0 m水位時減弱達到了9.6%，最大相對位移在高水位時減弱效果要比低水位更加明顯，在水位高度4.0 m時差值最大。

(3) SPH 法模擬各水位工況對最大相對位移、最大加速度和樓層加速度反應譜的結果表明：水箱水體的晃動對安全殼水箱以上部位的影響要大于對安全殼下部結構的影響；水位高度4.0 m 時，最大位移和最大加速度都出現了最小值，樓層加速度反應譜也在4.0 m 水位時最小，結果表明：這一水位高度既能滿足較高的儲水量，又兼具不錯的減震效果，是該型安全殼水箱較為合理的儲水水位。