U-50Zr螺旋十字燃料熱力耦合特性分析

蔡孟珂， 叢騰龍， 顧漢洋

(上海交通大學(xué) 機(jī)械與動(dòng)力工程學(xué)院， 上海 200240)

為了提升現(xiàn)有壓水堆的功率水平和安全裕量，研究者多從燃料材料和燃料幾何2個(gè)角度對(duì)傳統(tǒng)核燃料進(jìn)行改進(jìn)[1-5]。在材料方面，U-50Zr合金材料具有較高的熱導(dǎo)率、改進(jìn)的抗輻照性能、較好的包殼相容性以及在閉式循環(huán)中的優(yōu)勢(shì)，因此極具潛力[6]。在核燃料服役的過(guò)程中，氣體裂變產(chǎn)物的堆積使得燃料晶胞內(nèi)和晶界間形成微小氣泡或孔隙，進(jìn)而改變?nèi)剂系暮暧^熱力性質(zhì)。目前，針對(duì)傳統(tǒng)UO2陶瓷燃料的孔隙效應(yīng)的研究已較為完備；然而，針對(duì)U-50Zr合金材料的研究卻幾乎為空白。為了準(zhǔn)確地預(yù)測(cè)U-50Zr合金核燃料的堆內(nèi)輻照性能，需要對(duì)其材料孔隙效應(yīng)進(jìn)行研究。在燃料幾何方面，螺旋十字燃料(helical cruciform fuel，HCF)是一種結(jié)合了棒狀與板狀燃料元件優(yōu)點(diǎn)的革新型設(shè)計(jì)。HCF燃料的橫截面為十字型，其翼片沿軸向方向扭轉(zhuǎn)形成螺旋結(jié)構(gòu)，可提高燃料的換熱面積，并增強(qiáng)流體的攪渾換熱性能。同時(shí)，燃料的翼根區(qū)離燃料中心的徑向距離較短，減小芯體導(dǎo)熱路徑，降低燃料中心溫度。另外，在軸向每扭轉(zhuǎn)1/4周期平面內(nèi)，相鄰燃料間的翼尖區(qū)相互接觸，形成“自定位”結(jié)構(gòu)，故可省略定位格架，并在一定程度上彌補(bǔ)燃料格架帶來(lái)的壓降損失。

在新型核燃料的設(shè)計(jì)和優(yōu)化中，燃料的熱工水力特性以及熱力響應(yīng)反饋分析極為重要。目前，對(duì)于HCF的研究主要聚焦于熱工水力領(lǐng)域，主要包括HCF的壓降損失，攪渾特性以及傳熱能力。文獻(xiàn)[7-8]通過(guò)實(shí)驗(yàn)和數(shù)值模擬的方法，對(duì)HCF的壓降特性進(jìn)行了研究。文獻(xiàn)[9-14]通過(guò)視網(wǎng)探針和CFD模擬技術(shù)對(duì)HCF棒束的橫向攪渾性能進(jìn)行了研究，并預(yù)測(cè)了其在單相和兩相情況下的傳熱特性。李惠悅等[15]采用SSTk-ω湍流模型對(duì)不同螺距的星形螺旋燃料棒束通道內(nèi)流動(dòng)換熱特性進(jìn)行數(shù)值模擬研究。Shirvan[16]基于兩相歐拉方法研究了高壓下HCF的沸騰危機(jī)現(xiàn)象。除了上述熱工水力研究之外，燃料熱力性能分析也對(duì)核燃料的設(shè)計(jì)和優(yōu)化具有重要意義。然而，目前針對(duì)HCF的燃料性能分析的研究十分稀少，鄧陽(yáng)斌等[17]對(duì)單根HCF的熱力性能分析進(jìn)行了初步預(yù)測(cè)，但是，由于HCF復(fù)雜的三維幾何結(jié)構(gòu)，需要對(duì)多棒束下HCF的熱力響應(yīng)特性進(jìn)行更加深入和全面的研究，以捕捉包殼翼尖處的接觸應(yīng)力特性。

本文從材料和燃料分析2方面出發(fā)，對(duì)U-50Zr合金材料HCF進(jìn)行多物理場(chǎng)耦合分析。在材料方面，本文運(yùn)用分子動(dòng)力學(xué)軟件LAMMPS，建立了帶空隙U-50Zr合金材料的半經(jīng)驗(yàn)修正模型。在燃料分析方面，本文考慮了包括U-50Zr材料孔隙退化模型在內(nèi)的多物理場(chǎng)效應(yīng)，對(duì)較深輻照下HCF的熱力響應(yīng)特性進(jìn)行了耦合模擬和分析。

1 熱力耦合數(shù)學(xué)物理模型

本文從材料和燃料分析2方面出發(fā)，對(duì)U-50Zr合金材料HCF進(jìn)行多物理場(chǎng)耦合分析，如圖1所示。將分別給出U-50Zr合金彈性模量、熱導(dǎo)率、熱力耦合計(jì)算模型。

1.1 楊氏模量的計(jì)算

為了獲取不同溫度下的U-50Zr合金材料的彈性模量，常用的方法是進(jìn)行單軸拉伸實(shí)驗(yàn)?zāi)M。對(duì)大多數(shù)材料而言，楊氏模量可看作是各向同性的。然而，U-50Zr合金是由δ相UZ2組成的，該δ相UZ2是六方點(diǎn)陣排布的[18]，因此U-50Zr合金的彈性模量是各向異性的。為了獲取更準(zhǔn)確的彈性模量預(yù)測(cè)值，單軸拉伸實(shí)驗(yàn)需在相互正交的3個(gè)方向上進(jìn)行。在三維正交坐標(biāo)系下，應(yīng)力和應(yīng)變的關(guān)系為：

(1)

式中：下標(biāo)11、22、33、12、13和23分別表示x、y、z、xy、xz、yz方向；Cij為彈性模量張量。為了便于工程使用，往往將該6×6彈性模量矩陣轉(zhuǎn)化成等效楊氏模量，其轉(zhuǎn)化關(guān)系為：

(2)

(3)

(4)

(5)

(6)

(7)

式中：K是體積模量；G是剪切；E是楊氏模量。

圖1 基于U-50Zr合金材料下螺旋十字燃料的多物理場(chǎng)耦合Fig.1 The multi-physics coupling analysis for HCF based on U-50Zr alloy

1.2 熱導(dǎo)率的計(jì)算

U-50Zr合金的熱導(dǎo)率kt由聲子導(dǎo)熱率kp和電子導(dǎo)熱率ke2部分組成：

kt=kp+ke

(8)

本文采用非平衡分子動(dòng)力學(xué)方法預(yù)測(cè)材料的聲子導(dǎo)熱率，其原理由圖2所示。在長(zhǎng)平板的兩端施加較冷的溫度邊界，在平板的中段施加較熱的溫度邊界，通過(guò)測(cè)量平板材料的溫度梯度和熱流密度，進(jìn)而由傅里葉導(dǎo)熱定律求出聲子導(dǎo)熱率為：

(9)

式中：vh和vc是熱端和冷端的原子平均速度；t是導(dǎo)熱時(shí)間；S是平板端面截面積大??；dT/dx是溫度梯度。

圖2 非平衡分子動(dòng)力學(xué)原理求聲子熱導(dǎo)率示意Fig.2 The sketch of the NEMD method to obtain phonon conductivity

同時(shí)，電子導(dǎo)熱率跟電子阻抗值有關(guān)，二者的關(guān)系可通過(guò)Wiedemann-Franz定律描述：

(10)

式中：L是洛倫茲常數(shù)；ρe是電子阻抗。本文電子阻抗的計(jì)算借用了Zhou等[19]提出的經(jīng)驗(yàn)?zāi)Ｐ停?/p>

(11)

式中：B是電子-聲子散射時(shí)間；C是電子-電子散射時(shí)間；ρres是殘余電阻率；ρsat是飽和電阻率。

1.3 熱力耦合計(jì)算

熱力耦合計(jì)算由瞬態(tài)熱傳導(dǎo)方程和力學(xué)本構(gòu)方程組成。其中，熱傳導(dǎo)方程由時(shí)間項(xiàng)，擴(kuò)散項(xiàng)和源項(xiàng)構(gòu)成：

(12)

(13)

(14)

式中：σ是柯西應(yīng)力；Fv是單位體積力；E是彈性模量；v是泊松比；ε是彈性應(yīng)變。在非彈性變形階段，總應(yīng)變由彈性應(yīng)變、熱應(yīng)變、塑性應(yīng)變、腫脹應(yīng)變和蠕變應(yīng)變5部分組成。為了實(shí)現(xiàn)完全熱力耦合，導(dǎo)熱和變形方程同時(shí)求解，并同時(shí)考慮熱對(duì)力、力對(duì)熱的影響，故全局控制方程為：

(15)

式中：θ是節(jié)點(diǎn)溫度；u是節(jié)點(diǎn)位移；Kθθ是溫度剛度矩陣；Kuu是位移剛度矩陣；Kθu是位移對(duì)溫度影響的剛度矩陣；Kuθ是溫度對(duì)位移影響的剛度矩陣；Rθ和Ru分別是溫度和位移方程的殘差。

2 熱力耦合結(jié)果分析

2.1 材料孔隙效應(yīng)分析

本文預(yù)測(cè)了U-50Zr合金在溫度范圍為300～800 K的等效楊氏模量。結(jié)合裂變氣泡在燃料中的尺寸大小，孔隙尺寸半徑設(shè)置為0.5～1.2 nm，對(duì)應(yīng)的孔隙率大小為0.37%～7.99%，詳見(jiàn)表1。圖3給出了不同溫度、不同孔隙率下，楊氏模量的變化趨勢(shì)圖。由圖可見(jiàn)，隨著溫度和孔隙率的升高，U-50Zr合金的楊氏模量逐漸減小。

表1 楊氏模量計(jì)算中孔隙半徑和孔隙率的設(shè)置

圖3 不同溫度和孔隙率下U-50Zr合金的楊氏模量云圖Fig.3 The contour of Young’s modulus of U-50Zr at different temperatures and porosities

同理，本文預(yù)測(cè)了U-50Zr合金在溫度范圍為300～800 K的熱導(dǎo)率，孔隙尺寸和孔隙率的設(shè)置情況見(jiàn)表2。圖4給出了不同溫度、不同孔隙率下，熱導(dǎo)率的變化趨勢(shì)圖。由圖4可見(jiàn)，隨著孔隙率的增大，U-50Zr合金的熱導(dǎo)率逐漸降低；但是，隨著溫度的升高，U-50Zr合金的熱導(dǎo)率逐漸增高，這是因?yàn)楦邷靥岣吡穗娮訉?dǎo)熱率。

表2 熱導(dǎo)率計(jì)算中孔隙半徑和孔隙率的設(shè)置

圖4 不同溫度和孔隙率下U-50Zr合金的熱導(dǎo)率云圖Fig.4 The contour of thermal conductivity of U-50Zr at different temperatures and porosities

為了更加準(zhǔn)確地表述楊氏模量、熱導(dǎo)率和孔隙率的關(guān)系，本文建立了帶空隙U-50Zr合金材料的半經(jīng)驗(yàn)修正模型為：

Ep=(150.1-0.052 7T)(1-9.5P)

(16)

kp=(2.245 71-7.73×10-3T+

1.071 43×10-5T2)(1-2.15P)

(17)

該空隙效應(yīng)模型可初步預(yù)測(cè)U-50Zr合金的楊氏模量和熱導(dǎo)率，其適用于低空隙率(p<10%)和納米級(jí)氣泡(氣泡半徑為納米量級(jí))。

2.2 燃料熱力響應(yīng)分析

本研究基于有限元軟件Abaqus的Standard隱式求解器，綜合考慮基于U-50Zr的HCF多物理場(chǎng)現(xiàn)象：熱膨脹、輻照腫脹、蠕變、彈塑性變形等，并植入了式(16)和式(17)中的孔隙效應(yīng)模型，對(duì)深輻照下HCF的熱力響應(yīng)特性進(jìn)行分析。HCF在線(xiàn)功率水平18.2 kW/m下輻照了2.4 a。若以發(fā)生裂變的原子數(shù)占總原子數(shù)的百分比(per fissions of initial mental atoms, FIMA)為燃耗單位，其燃耗深度達(dá)到了14.1% FIMA。圖5給出了典型HCF燃料的幾何結(jié)構(gòu)。另外，由于HCF具有翼片結(jié)構(gòu)，因此翼根區(qū)域的導(dǎo)熱路徑小于翼尖區(qū)域，因此大部分熱量從翼根區(qū)域出流出，造成翼根處的熱流密度遠(yuǎn)大于翼尖處。

圖5 螺旋十字燃料3維幾何構(gòu)型Fig.5 The 3D geometry of helical cruciform fuel

在傳熱方面，圖6和圖7分別給出了未輻照情況下，HCF的溫度和熱流密度分布云圖。由于U-50Zr體系下的HCF具有較高的熱導(dǎo)率和換熱性能，因此HCF的整體溫度偏低，在軸向250 mm處的活性區(qū)頂端，最高溫度為674.4 K，遠(yuǎn)低于常規(guī)陶瓷棒狀燃料。

圖6 未輻照下螺旋十字燃料的溫度云圖Fig.6 Temperature contour of unirradiated HCF

圖7 未輻照下螺旋十字燃料的熱流密度云圖Fig.7 Heat flux contour of unirradiated HCF

圖8為包殼外壁面的溫度和熱流密度的周向分布圖。顯然，在45°、135°、225°、315°(即翼根)處的溫度和熱流密度均大于翼尖處。其中，溫度在翼根處達(dá)到最大值583 K，在翼尖處達(dá)到最小值564 K；熱流密度在翼根處達(dá)到最大值1 150 kW/m2,在翼尖處達(dá)到最小值281 kW/m2。圖9為燃料中心溫度隨燃耗變化的影響。隨著燃耗的加深，燃料產(chǎn)生中的裂變氣體不斷積累，進(jìn)而使燃料產(chǎn)生孔隙，減小燃料的熱導(dǎo)率。但是由圖9看出，燃料中心溫度隨燃耗的變化不明顯，從未經(jīng)輻照到輻照末期，溫度提升不超過(guò)0.3 K，這充分體現(xiàn)了U-50Zr體系下HCF的耐輻照穩(wěn)定性。

圖8 未輻照下螺旋十字燃料溫度和熱流密度周向分布圖Fig.8 Circumferential distribution of temperature and heat flux of unirradiated HCF

在受力方面，圖10 (a)給出了未輻照下，軸向扭轉(zhuǎn)90°平面內(nèi)HCF燃料的Mise應(yīng)力分布情況，應(yīng)力最大值位于產(chǎn)生接觸的包殼翼尖區(qū)域。由圖10 (b)所示，當(dāng)燃耗達(dá)到14.1% FIMA時(shí)，應(yīng)力最大值轉(zhuǎn)移到包殼翼尖區(qū)域兩側(cè)，且應(yīng)力集中現(xiàn)象較未輻照時(shí)有所緩解，應(yīng)力峰值也相應(yīng)降低。圖11給出了軸向扭轉(zhuǎn)90°平面內(nèi)HCF燃料翼尖區(qū)域應(yīng)力峰值隨燃耗的變化情況。在燃耗從0～1.1% FIMA內(nèi)，由于蠕變效應(yīng)的影響，應(yīng)力峰值從326.7 MPa下降到307.3 MPa。隨著燃耗進(jìn)一步的深入，輻照腫脹加大了燃料向外膨脹的趨勢(shì)，因此使得接觸處的應(yīng)力峰值重新開(kāi)始增大，當(dāng)燃耗達(dá)到3.3% FIMA 時(shí)，應(yīng)力峰值增大到322.1 MPa。此后，燃料應(yīng)力在蠕變的影響下開(kāi)始回落。最終，當(dāng)燃耗達(dá)到14.1% FIMA 時(shí)，翼尖處應(yīng)力峰值達(dá)到了313.3 MPa。在相同情況下，Geelhood等[20]提出的經(jīng)驗(yàn)?zāi)Ｐ陀?jì)算得到鋯合金包殼的極限斷裂強(qiáng)度值為500 MPa左右。因此，在較深燃耗下，U-50Zr體系下HCF的包殼應(yīng)力最大值低于材料失效強(qiáng)度，表明其擁有較好的力學(xué)安全性能。

圖9 螺旋十字燃料中心溫度隨燃耗變化Fig.9 Plot of central temperature of HCF

圖10 螺旋十字燃料的Mises應(yīng)力云圖Fig.10 Mises stress contour of HCF

圖11 螺旋十字燃料應(yīng)力峰值隨燃耗變化圖Fig.11 Plot of peak stress of HCF under different burnups

3 結(jié)論

1)在低孔隙率、低溫條件下，U-50Zr合金材料的楊氏模量隨孔隙率、溫度的升高而降低；U-50Zr合金材料的熱導(dǎo)率隨孔隙率的升高而降低，隨溫度的升高而升高。

2)螺旋十字燃料存在周向溫度、熱流密度分布的不均勻現(xiàn)象。

3)螺旋十字燃料的燃料溫度較低，在整個(gè)計(jì)算壽期內(nèi)，燃料中心溫度較為穩(wěn)定，且不超過(guò)675 K；另外，燃耗深度對(duì)溫度的影響不顯著。

4)U-50Zr材料的物性孔隙模型的準(zhǔn)確性仍需要未來(lái)的輻照實(shí)驗(yàn)數(shù)據(jù)的驗(yàn)證和支撐。