螺旋十字燃料棒束熱力耦合特性研究

叢騰龍 劉峪潔 郭 輝 肖 瑤 顧漢洋

1（上海交通大學(xué) 核科學(xué)與工程學(xué)院 上海 200240）

2（哈爾濱工程大學(xué) 核科學(xué)與技術(shù)學(xué)院 哈爾濱 150009）

燃料組件是反應(yīng)堆的核心釋熱部件，是決定反應(yīng)堆安全性和經(jīng)濟(jì)性的關(guān)鍵要素之一，也是反應(yīng)堆研發(fā)和設(shè)計(jì)的重點(diǎn)。螺旋十字燃料（Helical-Cruciform Fuel，HCF）是一種革新型的燃料組件，HCF元件具體結(jié)構(gòu)如圖1所示，十字形橫截面結(jié)構(gòu)沿軸向旋扭形成螺旋十字，通過緊密排布，可在高度方向上形成周期性自支撐結(jié)構(gòu)，無須定位格架。同時(shí)，十字截面提高了堆芯比換熱面積、降低導(dǎo)熱距離，且可以通過螺旋結(jié)構(gòu)增強(qiáng)通道間擾流交混，可在維持堆芯安全裕量的前提下提高堆芯功率密度，用于先進(jìn)緊湊型堆芯的設(shè)計(jì)。目前俄羅斯已經(jīng)將具有螺旋翼片的燃料組件用于高溫氣冷堆中［1］，美國(guó)也在積極推進(jìn)HCF組件的工程應(yīng)用，圖2為美國(guó)Lightbridge公司在美國(guó)能源部的資助下提出的正方形柵格HCF組件設(shè)計(jì)方案，與傳統(tǒng)17×17圓棒組件在尺寸上兼容，可以提供更高的單盒組件功率輸出［2］。

圖1 單根螺旋十字燃料及自支撐結(jié)構(gòu)示意圖(a) 單根燃料，(b) 自支撐結(jié)構(gòu)Fig.1 Sketch for the single fuel rod and the self-support structure (a) Single fuel rod, (b) Self-support structure

圖2 Lightbridge公司設(shè)計(jì)用于壓水堆的HCF組件截面圖Fig.2 Cross section of the HCF assembly for PWR designed by Lightbridge

現(xiàn)有針對(duì)HCF組件的研究多集中在流動(dòng)傳熱方面，包括單相阻力［3-5］、交混［3，6］、對(duì)流傳熱［5］、沸騰傳熱［6］和臨界熱流密度［4，6］，證明了HCF組件相比傳統(tǒng)帶格架圓棒組件在流阻、交混、傳熱和臨界熱流密度方面具有優(yōu)勢(shì)。然而，HCF組件依靠自身接觸支撐形成自定位結(jié)構(gòu)，在反應(yīng)堆運(yùn)行過程中，燃料相對(duì)冷態(tài)會(huì)出現(xiàn)明顯的熱膨脹，棒束間點(diǎn)接觸區(qū)域可能會(huì)出現(xiàn)應(yīng)力集中，局部過高的應(yīng)力會(huì)使包殼出現(xiàn)塑性變形乃至破裂，危害燃料組件的完整性，MIT的Deng等［7］對(duì)輻照條件下單根燃料的應(yīng)力應(yīng)變進(jìn)行了分析，通過引入基體和包殼在輻照作用下的相關(guān)熱、力學(xué)模型，得到金屬基體的HCF單棒使役性能，但由于僅考慮單棒，對(duì)于多棒之間接觸及其對(duì)應(yīng)力應(yīng)變影響規(guī)律的考慮仍不充分。鑒于此，本文建立典型HCF棒束單元的熱力耦合分析模型，研究了正常穩(wěn)瞬態(tài)和事故工況下HCF元件的應(yīng)力應(yīng)變分布，評(píng)價(jià)了HCF組件抵御事故的能力，為HCF組件設(shè)計(jì)提供了參考。

1 螺旋十字燃料熱力耦合數(shù)學(xué)物理模型

相對(duì)于冷態(tài)停堆工況，反應(yīng)堆熱態(tài)運(yùn)行階段HCF溫度會(huì)上升200~400 K，產(chǎn)生熱膨脹和熱應(yīng)力，HCF棒束翼片頂部的接觸位置出現(xiàn)應(yīng)力集中，在熱應(yīng)力和局部集中應(yīng)力共同作用下燃料包殼可能出現(xiàn)大幅度的彈性和塑性變形。此外，十字截面導(dǎo)致的燃料棒內(nèi)部非均勻徑向?qū)崾沟萌剂习魞?nèi)部溫度差異很大，影響材料的熱力學(xué)物性，進(jìn)一步復(fù)雜化HCF元件內(nèi)部應(yīng)力應(yīng)變場(chǎng)。本文分別建立燃料元件的溫度場(chǎng)和結(jié)構(gòu)場(chǎng)控制方程，給定熱對(duì)流和約束邊界條件，通過熱-力耦合求解熱力方程獲得燃料元件的傳熱和力學(xué)響應(yīng)。

1.1 基本控制方程

1.1.1 溫度場(chǎng)方程

螺旋十字燃料組件的溫度場(chǎng)可看作含有內(nèi)熱源的導(dǎo)熱問題，包殼和芯塊兩個(gè)固體區(qū)域的溫度分布可以由導(dǎo)熱微分方程確定：

穩(wěn)態(tài)工況下式（1）可以改寫為：

式中：ρ為密度，kg?m-3；c為比熱容，J?（kg·K）-1；T為溫度，K；τ為時(shí)間，s；λ為導(dǎo)熱系數(shù)，W?（m·K）-1；φ為單位時(shí)間內(nèi)單位體積中產(chǎn)生的熱量，W?m-3，對(duì)于包殼，該項(xiàng)為零。

1.1.2 結(jié)構(gòu)力學(xué)場(chǎng)控制方程

反應(yīng)堆運(yùn)行階段螺旋十字燃料組件的溫度會(huì)升高，由于燃料間的自支撐結(jié)構(gòu)限制了燃料的熱膨脹，燃料內(nèi)部產(chǎn)生熱應(yīng)力，均質(zhì)各向同性體滿足的基本方程為：

式中：σij為應(yīng)力張量，Pa；εij為應(yīng)變張量；ui為位移張量，m；ν為泊松比；E為彈性模量，Pa；δij為kronecker符號(hào)；α為熱膨脹系數(shù)，K-1；ΔT為與參考溫度的溫差，K；εpl為塑性應(yīng)變張量，在彈性階段為零。本文采用von Mises屈服準(zhǔn)則判斷材料是否進(jìn)入塑性階段，有效應(yīng)力σs可以表示為：

式中：σ1、σ2、σ3分別為第一、第二、第三主應(yīng)力。當(dāng)材料進(jìn)入塑性階段后，根據(jù)Levy-Mises增量理論，εpl可以寫成：

1.2 燃料的材料模型

由于U-50wt%Zr合金的熱導(dǎo)率高、裂變?cè)用芏却笄逸椪漳[脹率低，螺旋十字燃料采用U-50wt%Zr合金作為芯塊材料，同時(shí)采用Zr-4合金作為包殼材料。包殼和芯塊通過冶金結(jié)合的方式連接，中間無氣隙。

1.2.1 芯塊

在一些嚴(yán)重事故中芯塊中心溫度可能超過800 K，此時(shí)U-50wt%Zr合金會(huì)發(fā)生相變，導(dǎo)致材料的物性出現(xiàn)變化，但國(guó)內(nèi)外對(duì)相變后U-50wt%Zr合金的研究較少，且本文研究工況范圍內(nèi)芯塊溫度均低于相變點(diǎn)溫度（800 K），故不考慮相變對(duì)物性的影響。U-50wt%Zr合金的熱膨脹系數(shù)αf為［8］：

芯塊的密度［9］為 9.5×103kg?m-3，導(dǎo)熱系數(shù)和比熱分別依據(jù) Bauer［9］和 Fedorov［10］的研究結(jié)果計(jì)算。根據(jù)修正后的混合物定律［11］，芯塊的彈性模量Ef和泊松比νf可以表示為：

式中：EU為588 K下純鈾的彈性模量，EU=1.6×1011Pa；P為燃料的孔隙率，對(duì)于金屬芯體的燃料，在未經(jīng)輻照條件下孔隙率為零；WZ為Zr的質(zhì)量百分比，本文中WZ=0.5；TmU為純鈾的熔點(diǎn)，TmU=1 405 K；νU為588 K下純鈾的泊松比，νU=0.24。

1.2.2 包殼

當(dāng)300 K＜T＜1 800 K時(shí)，包殼導(dǎo)熱系數(shù)λc的表達(dá)式為［12］：

包 殼 的 密 度［9］為 6.55×103kg ?m-3，采 用MATPRO［12］物性手冊(cè)中的模型計(jì)算包殼的比熱，包殼的熱膨脹系數(shù)［13］與溫度無關(guān)，為5.58×10-6K-1，彈性模量Ec和泊松比νc采用Fisher模型［13］：

進(jìn)入塑性階段后包殼的應(yīng)力強(qiáng)化曲線可以擬合成［13］：

式中：σtrue和εtrue為真實(shí)應(yīng)力和真實(shí)應(yīng)變?yōu)檎鎸?shí)應(yīng)變率，如

1.3 有限元模型

圖2所示的燃料組件包括多根燃料棒，對(duì)整盒組件進(jìn)行力學(xué)分析所需的計(jì)算代價(jià)極大，由圣維南定律可知，組件外殼的約束僅影響邊緣燃料肋片頂部區(qū)域的受力，中心燃料的力學(xué)行為主要與相鄰燃料的約束有關(guān)。考慮到本文的目的是研究中心螺旋十字燃料棒在不同傳熱學(xué)邊界條件下的熱力響應(yīng)，因此，僅以包含中心目標(biāo)棒及周圍的3×3棒束單元作為計(jì)算域進(jìn)行熱力耦合分析。此外，考慮到3×3棒束中4根角棒與中心棒之間無直接力學(xué)作用，分析中忽略4根角棒，僅建立圖3所示的5根棒的計(jì)算域模型，單根燃料棒總長(zhǎng)為0.8 m，其中燃料段長(zhǎng)0.6 m，兩端的漸變段和固定段各為0.05 m，通過周圍4根棒為中間棒提供較為真實(shí)的約束條件。

圖3 計(jì)算域模型 (a) 燃料組件的3D視圖，(b) 扭轉(zhuǎn)角90°、180°、270°截面圖，(c) 單根燃料幾何尺寸Fig.3 Computational domain of the rod bundle (a) 3D view of the rod bundle, (b) Cross section at the 90°, 180°, 270° plane,(c) Dimension of a single rod

1.3.1 網(wǎng)格劃分

采用六面體網(wǎng)格對(duì)HCF組件進(jìn)行網(wǎng)格劃分。單根燃料在包殼外表面沿周向存在184個(gè)節(jié)點(diǎn)，其中肋片頂部區(qū)域的網(wǎng)格較密，可以提高接觸計(jì)算的精度。從燃料中心到肋片頂部徑向共有13個(gè)節(jié)點(diǎn)，在肋片根部燃料的熱流密度較高，溫度梯度較大，為準(zhǔn)確計(jì)算肋片根部的溫度場(chǎng)，從燃料中心到包殼外表面沿徑向共劃分18個(gè)節(jié)點(diǎn)。各處包殼中沿徑向均有4個(gè)節(jié)點(diǎn)。HCF組件燃料段的軸向網(wǎng)格節(jié)點(diǎn)數(shù)為180，單根HCF棒束網(wǎng)格總數(shù)為317 088，共計(jì)374 000個(gè)節(jié)點(diǎn)，具體如圖4所示。

圖4 單根燃料的網(wǎng)格劃分(a) 單根燃料網(wǎng)格，(b) 軸向截面網(wǎng)格Fig.4 Mesh partitioning of a single fuel(a) Mesh of the single rod, (b) Mesh of the axial section

1.3.2 對(duì)流邊界條件

由于本文重點(diǎn)研究燃料的力學(xué)行為，為簡(jiǎn)化計(jì)算忽略燃料周圍的流體區(qū)域，在燃料組件的包殼外表面上給定第三類邊界條件，即對(duì)流傳熱系數(shù)和主流溫度。第三類邊界條件可表示為如下形式：

式中：n為壁面外法線；Tf為主流溫度，K；Tw為壁面溫度，K；h為壁面對(duì)流傳熱系數(shù)，W?（m2·K）-1，模擬單相液體對(duì)流或兩相沸騰傳熱工況時(shí)，分別采用Dittus-Boelter公式［14］和Chen公式［15］計(jì)算換熱系數(shù)，穩(wěn)態(tài)工況中流動(dòng)傳熱邊界條件見表1，在兩相工況的計(jì)算達(dá)到穩(wěn)定后引入破口失水事故和反應(yīng)性引入事故以分析燃料在更高的運(yùn)行參數(shù)下抵御事故的能力，瞬態(tài)事故工況的事故序列見表2。

表1 穩(wěn)態(tài)工況中的傳熱邊界Table 1 Heat transfer boundary under steady conditions

表2 瞬態(tài)事故工況的事故序列Table 2 Time sequence of the events for transient accidents

1.3.3 結(jié)構(gòu)約束條件

對(duì)HCF組件進(jìn)行力學(xué)模擬的外部約束條件如下：

1）燃料組件的四周設(shè)置一個(gè)無應(yīng)變的剛性外殼（圖3中6處），外殼的內(nèi)壁冷態(tài)時(shí)與邊緣棒發(fā)生點(diǎn)接觸，在接觸面的法向邊緣棒的位移等于零，邊緣棒與外殼和中心棒之間均存在摩擦，摩擦系數(shù)為0.1；

2）對(duì)燃料底部的固定段及圓形端面（圖3中5處）施加徑向、周向和軸向的固定約束，使該固定段的位移為零，模擬下管座的支撐和固定作用；

3）對(duì)燃料頂部的固定段（圖3中1處）施加徑向、周向的固定約束，在圓形端面設(shè)置軸向的彈簧約束，彈簧的剛度k為1.45×1010N?m-1，使該固定段在徑向和周向的位移等于零，軸向位移x受到彈簧反力F限制，具體方程為F+kx=0，模擬上管座的定位和壓緊作用。

1.4 耦合及數(shù)值求解方法

本文使用ANSYS Workbench進(jìn)行熱力耦合分析，由于燃料棒變形幅度較小，不考慮燃料棒變形對(duì)其內(nèi)部導(dǎo)熱的影響。先通過熱單元steady-state thermal（transient thermal）求解燃料的穩(wěn)態(tài)（瞬態(tài)）溫度場(chǎng)，再將數(shù)據(jù)單向傳遞至結(jié)構(gòu)單元static structural進(jìn)行穩(wěn)態(tài)（瞬態(tài)）力學(xué)分析。采用Newton-Raphson迭代法求解燃料在接觸區(qū)域和進(jìn)入塑性階段后的非線性方程，位移的容差因子為0.5%。

2 結(jié)果分析

本文對(duì)圖3所示的計(jì)算域開展數(shù)值模擬研究，獲得中心棒及周圍四根棒的熱力響應(yīng)，并重點(diǎn)對(duì)中心棒的熱力行為展開分析。螺旋十字燃料組件沿軸線發(fā)生扭轉(zhuǎn)，扭轉(zhuǎn)的角度與軸向高度成正比，為表明燃料組件的力學(xué)行為在軸向的周期性變化規(guī)律并對(duì)軸向高度進(jìn)行無量綱化，可定義扭轉(zhuǎn)角δ（rad）代替軸向高度z（0≤z≤0.6 m）表示燃料段的軸向位置，δ的表達(dá)式如下：

相鄰燃料在扭轉(zhuǎn)角為0°、90°、180°、270°和360°位置發(fā)生接觸。

2.1 單相液體對(duì)流傳熱工況

反應(yīng)堆正常運(yùn)行時(shí)，絕大部分冷卻劑處于單相液態(tài)，因此本文首先開展單相液體對(duì)流邊界條件下的螺旋十字燃料組件熱力行為研究，冷卻劑入口溫度為553.85 K，入口速度為1 m?s-1，燃料棒的體積釋熱率為740 MW?m-3，在冷卻劑入口溫度、速度和棒功率的基礎(chǔ)上，通過熱平衡計(jì)算冷卻劑溫度。利用Dittus-Bolter公式基于通道平均熱工參數(shù)計(jì)算第三類邊界條件所需的對(duì)流換熱系數(shù)，換熱系數(shù)沿棒周向?yàn)槌?shù)。通過計(jì)算獲得中心棒的von Mises應(yīng)力分布如圖5所示，可以看出相鄰燃料發(fā)生接觸的自支撐平面上肋片頂部區(qū)域的von Mises應(yīng)力較高，相對(duì)其他區(qū)域更可能發(fā)生塑性變形甚至破裂，因此本文主要對(duì)接觸平面上燃料的受力和變形進(jìn)行分析，其中扭轉(zhuǎn)角180°平面上的中心棒溫度、von Mises應(yīng)力和塑性應(yīng)變分布（同高度截面冷卻劑平均溫度為559.4 K）如圖6所示。燃料溫度最大值677 K處于芯塊中心處，周圍的等溫線形狀為圓形。而靠近包殼的燃料溫度由于受冷卻劑對(duì)流傳熱影響更大，等溫線形狀與燃料的十字形輪廓相似。等溫線形狀的變化使得燃料內(nèi)部徑向的溫度梯度沿周向出現(xiàn)較大差異，因溫度梯度產(chǎn)生的熱應(yīng)力在周向上也呈現(xiàn)不均勻分布（如圖6（b）所示）。中心棒與邊緣棒接觸的區(qū)域溫度梯度較小，但由于接觸約束的存在，受應(yīng)力集中作用，該處von Mises應(yīng)力出現(xiàn)極大值，超過了屈服極限，發(fā)生明顯的塑性變形（圖6（c））。從圖6（d）和（e）可以看出，燃料芯塊和包殼中彈性應(yīng)變和熱膨脹應(yīng)變的最大值僅為塑性應(yīng)變最大值的十分之一，因此本文僅對(duì)塑性應(yīng)變進(jìn)行分析。

圖5 中心棒von Mises應(yīng)力分布云圖Fig.5 Distribution of von Mises stress on the center rod

圖6 180°平面中心棒的溫度、von Mises應(yīng)力和應(yīng)變分布(a) 溫度分布，(b) von Mises應(yīng)力分布，(c) 塑性應(yīng)變分布，(d) 彈性應(yīng)變分布，(e) 熱膨脹應(yīng)變分布Fig.6 Temperature, von Misesstress and strain at the 180° plane (a) Temperature distribution, (b) von Mises stress distribution,(c) Plastic strain distribution, (d) Elastic strain distribution, (e) Thermal strain distribution

為進(jìn)一步量化相鄰棒約束作用和燃料內(nèi)部溫度梯度對(duì)中心棒應(yīng)力分量的影響，圖7給出了中心棒扭轉(zhuǎn)180°平面包殼外表面處三個(gè)應(yīng)力分量的分布曲線（平面內(nèi)0°位置見圖3（b））。圖中0°、90°、180°和270°的位置區(qū)域（以下稱翼尖處）受接觸約束影響更大，由于接觸約束主要限制燃料的徑向熱膨脹，對(duì)周向和軸向的熱膨脹影響較小，因此翼尖處的徑向應(yīng)力為-700 MPa左右的壓應(yīng)力，周向應(yīng)力和軸向應(yīng)力僅分別為徑向應(yīng)力的55%和70%。而45°、135°、225°和315°的位置區(qū)域（以下稱翼根處）離約束位置較遠(yuǎn)，受接觸約束影響較小，從圖6（a）可以看出，翼根處的溫度梯度較大，說明該處的熱應(yīng)力與溫度梯度相關(guān)，溫度梯度使得包殼外側(cè)的溫度較低，熱膨脹量比內(nèi)側(cè)小，限制了內(nèi)側(cè)的熱膨脹，導(dǎo)致包殼外側(cè)在周向和軸向受到反作用的拉應(yīng)力，數(shù)值分別可以達(dá)到250 MPa和130 MPa，由于文章中忽略了流體區(qū)域，不考慮流體壓力對(duì)燃料包殼的影響，徑向應(yīng)力接近于零。

圖7 180°平面包殼外表面處的應(yīng)力分量分布（加粗刻度線為零應(yīng)力線）Fig.7 Distributions of stress components on the outer surface of cladding at the 180° plane(the bold scale line represents the zero stress scale)

螺旋十字燃料棒在扭轉(zhuǎn)角為90°、180°和270°平面發(fā)生約束定位，不同扭轉(zhuǎn)角平面冷卻劑的平均溫度不同，圖8展示了不同扭轉(zhuǎn)角平面包殼外表面處von Mises應(yīng)力和塑性應(yīng)變的分布情況。可以看出，三個(gè)扭轉(zhuǎn)角平面包殼外表面處的von Mises應(yīng)力和塑性應(yīng)變曲線均以90°為周期變化，von Mises應(yīng)力在翼尖和翼根處較大，分別為300 MPa和216 MPa左右，塑性應(yīng)變?cè)谝砑馓庍_(dá)到最大值0.04，說明發(fā)生約束定位的平面包殼外表面處的von Mises應(yīng)力和塑性應(yīng)變的變化趨勢(shì)基本相同。但從局部放大圖中可以看出，在扭轉(zhuǎn)角為90°和270°平面上，von Mises應(yīng)力和塑性應(yīng)變曲線分別向逆時(shí)針方向和順時(shí)針方向偏轉(zhuǎn)，這是燃料軸向沿螺旋線向兩個(gè)固定端膨脹的結(jié)果。由于燃料頂部固定端的彈簧剛度較大，燃料所受約束在軸向幾乎關(guān)于扭轉(zhuǎn)角為180°平面對(duì)稱，因此這一平面上偏轉(zhuǎn)現(xiàn)象不明顯。此外區(qū)別于點(diǎn)接觸時(shí)的應(yīng)力應(yīng)變尖峰，單相工況下翼尖處出現(xiàn)了一段跨度約2°，峰值為300 MPa的von Mises應(yīng)力平臺(tái)和數(shù)值為0.02~0.04的塑性應(yīng)變平臺(tái)，說明翼尖處已經(jīng)從冷態(tài)時(shí)的點(diǎn)接觸轉(zhuǎn)變?yōu)槊娼佑|。

圖8 不同扭轉(zhuǎn)角平面包殼外表面處的von Mises應(yīng)力和塑性應(yīng)變分布(a) von Mises應(yīng)力（加粗刻度線為零應(yīng)力線），(b) 塑性應(yīng)變Fig.8 Distributions of von Mises stress and plastic strain on the cladding outer surface with different helical angles(a) von Misesstress (the bold scale line means the zero stress scale), (b) Plastic strain

2.2 沸騰傳熱工況

針對(duì)堆芯過熱通道中冷卻劑可能出現(xiàn)的流動(dòng)沸騰，本文研究了冷卻劑處于15.5 MPa下的飽和溫度617.9 K時(shí)體積釋熱率為740 MW?m-3的棒束的熱力響應(yīng)，由于計(jì)算區(qū)域加熱長(zhǎng)度較短，傳熱量較少，并且本文重點(diǎn)關(guān)注棒中間遠(yuǎn)離兩端影響的高度，即扭轉(zhuǎn)角為180°的位置，因此假設(shè)冷卻劑的質(zhì)量含氣率保持不變，為0.1。圖9給出了包殼內(nèi)外表面的溫度分布，可以看出由于兩相工況中冷卻劑的溫度較高，翼尖處和翼根處的溫度分別比單相工況中高54 K和38 K。但棒束的線功率密度保持不變，單相和兩相工況中翼尖處和翼根處包殼內(nèi)外表面的溫差相同，分別為5 K和20 K左右。

圖9 不同工況中180°平面包殼內(nèi)外表面的溫度分布Fig.9 Temperature distribution on the cladding inner and outer surfaces at 180° plane under different operating conditions

為研究溫度的差異對(duì)包殼力學(xué)行為的影響，圖10比較了單相和兩相工況扭轉(zhuǎn)180°平面包殼外表面處的von Mises應(yīng)力和塑性應(yīng)變分布。由于翼根處的熱應(yīng)力與包殼內(nèi)外表面的溫差有關(guān)，而兩種工況中溫差基本相同，翼根處的von Mises應(yīng)力在兩種工況中同為215 MPa。從圖10（b）可以看出，兩相工況翼尖處的塑性應(yīng)變較單相工況偏大15%~25%，說明翼尖處因熱膨脹產(chǎn)生的塑性應(yīng)變隨包殼溫度的上升而增加。但圖10（a）中兩相工況翼尖處的von Mises應(yīng)力卻較單相工況偏低6%~10%，這是由于溫度影響了包殼的力學(xué)性能。圖11展示了不同溫度下包殼在塑性階段的應(yīng)力強(qiáng)化曲線，可以看出溫度越高，包殼的屈服極限越低，說明相較于單相工況，兩相工況在較低的應(yīng)力水平下進(jìn)入塑性階段，塑性階段應(yīng)力隨應(yīng)變的增加速度較彈性階段放緩；圖11中曲線的斜率還隨著溫度的上升而減小，即同處于塑性階段時(shí)兩相工況中包殼的應(yīng)力隨塑性應(yīng)變的增加速度比單相工況中慢。這兩點(diǎn)原因共同導(dǎo)致兩相工況中翼尖處von Mises應(yīng)力比單相工況小。

圖10 不同工況中180°平面包殼外表面的von Mises應(yīng)力和塑性應(yīng)變分布(a) von Mises應(yīng)力（加粗刻度線為零應(yīng)力線），(b) 塑性應(yīng)變Fig.10 Distributions of von Mises stress and plastic strain on the cladding inner and outer surfaces at 180° plane under different operating conditions (a) von Mises stress (the bold scale line means the zero stress scale), (b) Plastic strain

圖11 包殼材料在不同溫度下應(yīng)力強(qiáng)化曲線Fig.11 Stress intensification curves for cladding material under different temperature

2.3 反應(yīng)性引入事故

本文研究的反應(yīng)性引入事故為正常運(yùn)行時(shí)控制棒突然提出堆芯，導(dǎo)致堆芯出現(xiàn)較大的正反應(yīng)性的事故，即彈棒事故。由圖10可知，單相和兩相工況中扭轉(zhuǎn)角180°平面包殼的應(yīng)力應(yīng)變最大值均出現(xiàn)在翼尖處和翼根處，出于安全分析的保守性考慮，本文對(duì)事故過程中包殼翼根處和翼尖處的力學(xué)行為進(jìn)行研究并于包殼的內(nèi)外表面設(shè)置了4個(gè)檢測(cè)點(diǎn)，如圖12（a）所示。反應(yīng)性引入事故于1 s時(shí)在處于額定功率正常運(yùn)行的反應(yīng)堆中發(fā)生，1.1 s時(shí)控制棒完全彈出堆芯，導(dǎo)致燃料功率在1.24 s達(dá)到峰值，為額定功率的1.83倍，包殼的溫度也隨之上升，監(jiān)測(cè)點(diǎn)中4處的溫度最高，最大值出現(xiàn)在1.6 s，為662 K，可以將事故以該時(shí)刻為界分為前期和后期兩個(gè)階段。由于反應(yīng)堆存在固有安全性，1.24 s后燃料功率不會(huì)繼續(xù)增加。與此同時(shí)收到1.1 s發(fā)出的停堆信號(hào)后，控制棒在2 s時(shí)開始下落并在5 s時(shí)完全插入堆芯，燃料功率迅速下降，5.5 s時(shí)僅為額定功率的20%，包殼溫度也降低至630 K以下，圖12（b）為反應(yīng)性引入事故中燃料功率［16］與各監(jiān)測(cè)點(diǎn)溫度的變化情況。

圖12 反應(yīng)性引入事故中燃料的功率和各監(jiān)測(cè)點(diǎn)溫度的變化 (a) 監(jiān)測(cè)點(diǎn)位置，(b) 功率和溫度的變化Fig.12 Variations of fuel power and temperature at different monitoring points(a) Locations of the monitoring points, (b) Power and temperature variations

反應(yīng)性引入事故中溫度的起伏變化會(huì)極大影響包殼的受力和變形，圖13展示了翼尖處包殼內(nèi)外表面的應(yīng)力隨時(shí)間的變化情況?？梢园l(fā)現(xiàn)，由于接觸約束限制了包殼徑向的熱膨脹，翼尖處徑向應(yīng)力的變化幅度明顯高于周向應(yīng)力和軸向應(yīng)力，其中包殼內(nèi)表面的徑向應(yīng)力在1 s和7.5 s時(shí)的差值分別是周向應(yīng)力和軸向應(yīng)力的數(shù)百倍和4.6倍，包殼外表面則是3倍和1.7倍，說明翼尖處包殼的受力變化主要由徑向應(yīng)力引起。從von Mises應(yīng)力的變化趨勢(shì)中可以看出，2 s前翼尖處包殼內(nèi)外表面的von Mises應(yīng)力均基本保持不變，2 s后包殼內(nèi)表面的von Mises應(yīng)力始終隨著燃料溫度的降低而降低，7.5 s時(shí)不足100 MPa；而包殼外表面的von Mises應(yīng)力則是先下降再上升，最小值52 MPa出現(xiàn)在2.36 s，7.5 s時(shí)von Mises應(yīng)力上升至超過該溫度下的屈服極限281 MPa。出現(xiàn)這一差異主要是由于包殼經(jīng)歷的力學(xué)過程不同：在反應(yīng)性引入事故期間，燃料溫度先上升后下降，包殼從受熱膨脹變?yōu)槔鋮s收縮，發(fā)生先加載后卸載的過程，但由于塑性變形不會(huì)在卸載過程中減小，燃料溫度下降后塑性應(yīng)變無法恢復(fù)到事故前的水平，燃料冷卻收縮過程中翼尖處包殼外表面會(huì)出現(xiàn)反向加載的現(xiàn)象，徑向應(yīng)力在2 s后先逐漸降低至零然后反向增大，如圖13（b）所示。而包殼內(nèi)表面則始終處于卸載狀態(tài)，徑向應(yīng)力在2 s后逐漸降低，如圖13（a）所示。因此當(dāng)燃料溫度降低時(shí)包殼內(nèi)外表面的von Mises應(yīng)力呈現(xiàn)不同的變化趨勢(shì)。

圖13 反應(yīng)性引入事故中翼尖處包殼內(nèi)外表面應(yīng)力的變化 (a) 包殼內(nèi)表面，(b) 包殼外表面Fig.13 Stress at the cladding inner and outer surfaces of blade tip under reactivity insertion accident(a) Cladding inner surface, (b) Cladding outer surface

圖14進(jìn)一步給出了翼尖處包殼內(nèi)外表面塑性應(yīng)變隨時(shí)間的變化規(guī)律。事故前期翼尖處包殼內(nèi)表面的塑性應(yīng)變隨著溫度的增加上升了15%左右；包殼外表面的塑性應(yīng)變由于接觸約束限制基本保持不變。事故后期包殼內(nèi)表面的塑性應(yīng)變保持不變，為0.036 9，說明包殼內(nèi)表面僅發(fā)生卸載過程；而3 s后包殼外表面的塑性應(yīng)變則從0.033下降到0.027，說明3 s后包殼外表面處發(fā)生反向加載過程并產(chǎn)生反向的塑性變形，抵消了一部分事故前產(chǎn)生的塑性應(yīng)變，使得7.5 s時(shí)的塑性應(yīng)變相對(duì)于事故前降低。

圖14 反應(yīng)性引入事故中翼尖處包殼內(nèi)外表面塑性應(yīng)變的變化 (a) 包殼內(nèi)表面，(b) 包殼外表面Fig.14 Plastic strain at the cladding inner and outer surfaces of blade tip under reactivity insertion accident(a) Cladding inner surface, (b) Cladding outer surface

與接觸約束影響翼尖處的徑向應(yīng)力不同，翼根處的受力與包殼內(nèi)外表面的溫差有關(guān)，圖15展示了反應(yīng)性引入事故中翼根處包殼內(nèi)外表面應(yīng)力的變化情況。由于翼根處包殼內(nèi)外表面的溫差從1 s時(shí)的21 K降低到7.5 s時(shí)的3.7 K，受溫度梯度影響的周向應(yīng)力和軸向應(yīng)力也出現(xiàn)大幅下降，數(shù)值分別可以達(dá)到100~170 MPa和60~80 MPa，而徑向應(yīng)力的變化幅度則小于10 MPa，即徑向應(yīng)力的變化幅度遠(yuǎn)小于周向應(yīng)力和軸向應(yīng)力的變化幅度，von Mises應(yīng)力的變化主要由周向應(yīng)力和軸向應(yīng)力的變化引起。從圖15還可以看出，翼根處包殼內(nèi)外表面的周向應(yīng)力與軸向應(yīng)力2 s前上升幅度不超過20 MPa，2 s后均呈現(xiàn)下降的趨勢(shì)，沒有出現(xiàn)方向變化，并且從圖16可以看出，事故后期包殼翼根處的塑性應(yīng)變保持不變，說明翼根處在事故后期僅發(fā)生卸載過程，沒有出現(xiàn)反向加載。雖然翼根處的塑性應(yīng)變隨著溫度的增加在1.6 s時(shí)達(dá)到事故前的數(shù)倍，但由于翼根處的塑性應(yīng)變僅為翼尖處的1/10，因此從防止包殼破裂的角度考慮更需要關(guān)注翼尖處塑性應(yīng)變的變化。

圖15 反應(yīng)性引入事故中翼根處包殼內(nèi)外表面應(yīng)力的變化 (a) 包殼內(nèi)表面，(b) 包殼外表面Fig.15 Stress at the cladding inner and outer surfaces of blade elbow under reactivity insertion accident(a) Cladding inner surface, (b) Cladding outer surface

圖16 反應(yīng)性引入事故中翼根處包殼內(nèi)外表面塑性應(yīng)變的變化 (a) 包殼內(nèi)表面，(b) 包殼外表面Fig.16 Temperature at the cladding inner and outer surfaces of blade elbow under reactivity insertion accident(a) Cladding inner surface, (b) Cladding outer surface

2.4 一回路壓力邊界破口失水事故

在處于滿功率正常運(yùn)行的反應(yīng)堆中1 s時(shí)引入一回路壓力邊界破口失水事故，事故發(fā)生后冷卻劑從破口中大量流出，堆芯壓力迅速下降，1.6 s時(shí)降至12.6 MPa并觸發(fā)低壓停堆信號(hào)，忽略停堆響應(yīng)時(shí)間和控制棒下落時(shí)間，不考慮冷卻劑空泡效應(yīng)的影響，反應(yīng)堆功率在1.6 s前保持不變并在1.7 s降低至事故前的十分之一，破口失水事故中堆芯壓力［17］和燃料功率［18］的變化曲線如圖17所示。事故中假設(shè)一回路冷卻劑始終處于飽和沸騰階段，飽和溫度隨堆芯壓力降低而降低。

圖17 破口失水事故中堆芯壓力和功率的變化Fig.17 Variations of core pressure and power under the LOCA condition

為分析一回路壓力邊界破口失水事故中包殼翼尖處和翼根處的溫度、變形和受力的變化情況，仍選取§3.3中的4個(gè)監(jiān)測(cè)點(diǎn)作為研究對(duì)象（圖12（a）），圖18和圖19分別給出了破口失水事故中翼尖處各監(jiān)測(cè)點(diǎn)應(yīng)力和塑性應(yīng)變的變化曲線?？梢钥闯?，從事故發(fā)生到停堆包殼內(nèi)外表面的溫差基本不變，但由于冷卻劑的溫度隨堆芯壓力下降而降低且燃料棒功率保持不變，包殼內(nèi)外表面的溫度下降了10 K，結(jié)合圖19中停堆前包殼內(nèi)外表面的塑性應(yīng)變保持不變，說明停堆前翼尖處發(fā)生卸載過程，包殼內(nèi)外表面的von Mises應(yīng)力分別降低37 MPa和20 MPa。停堆后0.2 s內(nèi)，包殼內(nèi)外表面溫度下降的平均值同為10 K，但溫差則由于堆芯功率的劇烈下降而降低至停堆前的40%，包殼內(nèi)外表面von Mises應(yīng)力的變化幅度則可以達(dá)到停堆前的3.5倍和10倍，說明翼尖處包殼的受力不僅與限制熱膨脹的接觸約束有關(guān)，同時(shí)也受包殼內(nèi)外表面溫差的影響。從圖18中還可以看出，停堆后短時(shí)間內(nèi)包殼外表面的von Mises應(yīng)力就已超過屈服極限，而包殼內(nèi)表面的von Mises應(yīng)力在9.5 s時(shí)仍小于200 MPa，圖19中包殼外表面的塑性應(yīng)變由0.034 2下降至0.018 4而包殼內(nèi)表面的塑性應(yīng)變保持0.034不變，說明包殼外表面先發(fā)生卸載后出現(xiàn)反向加載而包殼內(nèi)表面僅存在卸載過程。

圖18 破口失水事故中翼尖處包殼內(nèi)外表面應(yīng)力的變化 (a) 包殼內(nèi)表面，(b) 包殼外表面Fig.18 Stress at the cladding inner and outer surfaces of blade tip under the LOCA condition(a) Cladding inner surface, (b) Cladding outer surface

圖19 破口失水事故中翼尖處包殼內(nèi)外表面塑性應(yīng)變的變化 (a) 包殼內(nèi)表面，(b) 包殼外表面Fig.19 Plastic strain at the cladding inner and outer surfaces of blade tip under the LOCA condition(a) Cladding inner surface, (b) Cladding outer surface

為研究一回路壓力邊界破口失水事故中包殼翼根處應(yīng)力和塑性應(yīng)變的變化情況，本文繪制了翼根處各監(jiān)測(cè)點(diǎn)應(yīng)力和塑性應(yīng)變隨溫度和包殼內(nèi)外表面溫差的變化曲線如圖20和圖21所示。從圖20可以看出，停堆前翼根處包殼內(nèi)外表面的溫度下降了6 K而溫差基本不變，von Mises應(yīng)力也基本不變。停堆后0.05 s內(nèi)，翼根處的溫度同樣下降了6 K左右，包殼內(nèi)外表面的溫差降低到停堆前的80%，von Mises應(yīng)力出現(xiàn)下降，說明翼根處離接觸位置較遠(yuǎn)，包殼整體溫度的下降不會(huì)使該位置發(fā)生卸載，應(yīng)力受包殼內(nèi)外表面溫差影響更大。停堆前由于溫差存在0.2 K的小幅度上升，圖21中可以觀察到包殼內(nèi)外表面的塑性應(yīng)變?cè)黾恿?%~3%，但仍遠(yuǎn)小于翼尖處。停堆后包殼的von Mises應(yīng)力不斷下降且未超過屈服極限，翼根處僅發(fā)生卸載過程，塑性應(yīng)變保持不變。

圖20 破口失水事故中翼根處包殼內(nèi)外表面應(yīng)力的變化 (a) 包殼內(nèi)表面，(b) 包殼外表面Fig.20 Stress at the cladding inner and outer surfaces of blade elbow under the LOCA condition(a) Cladding inner surface, (b) Cladding outer surface

圖21 破口失水事故中翼根處包殼內(nèi)外表面塑性應(yīng)變的變化 (a) 包殼內(nèi)表面，(b) 包殼外表面Fig.21 Plastic strain at the cladding inner and outer surfaces of blade elbow under the LOCA condition(a) Cladding inner surface, (b) Cladding outer surface

3 結(jié)語

本文對(duì)典型螺旋十字燃料HCF棒束單元進(jìn)行熱力耦合特性研究，模擬了包括單相對(duì)流、飽和沸騰、反應(yīng)性引入事故和失水事故在內(nèi)的多種工況，重點(diǎn)分析中心棒在自定位平面上的傳熱狀況和力學(xué)行為，評(píng)估了HCF組件的安全性能。根據(jù)計(jì)算結(jié)果，可以得到如下結(jié)論：

1） 所有工況下HCF組件中心棒包殼外表面處的von Mises應(yīng)力和塑性應(yīng)變最大值總是出現(xiàn)在相鄰燃料的接觸區(qū)域，其次是燃料翼片的翼根處。不同軸向位置發(fā)生約束定位的平面上包殼外表面在周向的受力與變形的變化趨勢(shì)相同，但相位存在偏差。翼尖處的應(yīng)力主要受接觸約束影響，同時(shí)也和包殼內(nèi)外表面溫差有關(guān)；翼根處的應(yīng)力則與包殼內(nèi)外表面的溫差有關(guān)。

2） 反應(yīng)堆正常運(yùn)行冷卻劑處于單相液態(tài)時(shí)燃料包殼翼尖處會(huì)出現(xiàn)0.02~0.04的塑性應(yīng)變平臺(tái)和峰值為300 MPa的von Mises應(yīng)力平臺(tái)，說明翼尖處已經(jīng)從點(diǎn)接觸變?yōu)槊娼佑|。

3） 相較于單相液體冷卻工況，飽和沸騰工況下包殼的平均溫度高35~60 K，但包殼內(nèi)外表面的溫差相同，導(dǎo)致翼根處von Mises應(yīng)力相等，翼尖處塑性應(yīng)變偏大15%~25%，但由于溫度越高包殼的屈服極限越低、塑性階段應(yīng)力隨溫度的升高對(duì)塑性應(yīng)變的增加越不敏感，翼尖處von Mises應(yīng)力反而偏小6%~10%。

4） 反應(yīng)性引入和失水事故工況下，包殼的von Mises應(yīng)力和塑性應(yīng)變分別小于350 MPa和0.04，且包殼外表面溫度低于鋯水反應(yīng)溫度；停堆后，翼尖處包殼外表面發(fā)生反向加載并再次進(jìn)入塑性階段，導(dǎo)致該處塑性應(yīng)變降低，von Mises應(yīng)力再次上升并超過屈服極限，需要注意這一反向應(yīng)力對(duì)燃料的二次損傷。

本文對(duì)HCF組件熱力耦合特性進(jìn)行了研究，得到了燃料組件的力學(xué)響應(yīng)，但由于HCF組件的熱流密度沿周向并非均勻分布，冷卻劑在同一軸向平面上的溫度也不盡相同，因此在熱力耦合的基礎(chǔ)上還需要對(duì)冷卻劑的流動(dòng)傳熱進(jìn)行模擬。此外，本文的熱力分析中僅考慮了壽期初的燃料組件，隨著反應(yīng)堆的運(yùn)行，燃料燃耗深度的增加會(huì)改變芯塊和包殼材料的物性，輻照對(duì)物性的改變以及輻照腫脹現(xiàn)象本文也未分析，可以增加冷卻劑計(jì)算模塊、建立燃料的輻照和燃耗模型使HCF組件的力學(xué)計(jì)算更加準(zhǔn)確。

作者貢獻(xiàn)聲明叢騰龍負(fù)責(zé)立題和設(shè)計(jì)計(jì)算工況，數(shù)據(jù)分析，論文修改；劉峪潔負(fù)責(zé)數(shù)值模擬，數(shù)據(jù)分析，論文初稿；郭輝數(shù)值模擬；肖瑤技術(shù)支持；顧漢洋獲取研究經(jīng)費(fèi)。