快堆帶繞絲棒束組件低雷諾數(shù)下的水力特性分析

程道喜，齊曉光，杜開文，翟偉明

(中國原子能科學研究院 反應堆工程技術研究部，北京 102413)

國際上的快堆組件大多采用帶繞絲棒束的形式，燃料棒可布置得很緊湊?？於讶剂辖M件棒束采用帶繞絲的結構定位的同時，可加強各子通道間冷卻劑的交混，強化燃料棒與冷卻劑的換熱，但繞絲的存在會增大棒束段的阻力。各國如美國、法國、日本、意大利等均針對燃料組件的水力特性(主要是棒束段的水力特性)進行了大量的實驗研究，并提出了一系列帶繞絲棒束組件的阻力壓降計算模型。Engel等[1]、Cheng和Todreas[2]利用實驗總結了壓降經(jīng)驗關系式。此外劉一哲等[3]在CRT模型和Engel等研究的基礎上，針對中國實驗快堆(CEFR)燃料組件的具體情況提出了ICRT壓降關系式，用以計算冷卻劑在湍流區(qū)、過渡流區(qū)和層流區(qū)的棒束壓降。

近年來，日益強大的計算能力使計算流體力學(CFD)手段得以應用到堆芯組件流動及傳熱分析的相關領域，如壓水堆組件格架熱工水力特性研究[4]。由于快堆燃料組件結構復雜，對于快堆燃料組件內(nèi)部流場及溫度場的CFD研究起步較晚。Ahmad等[5-14]利用CFD方法從不同角度進行了快堆帶繞絲棒束組件的研究。綜合各項研究結果，近些年國外學者利用CFD程序對快堆組件帶繞絲棒束的三維流動特性以及傳熱特性進行了大量的研究，但研究大多集中在湍流區(qū)，而對于處于過渡區(qū)甚至是層流區(qū)的研究很少。

層流及過渡流動這種低雷諾數(shù)(Re)流動工況對于快堆的安全運行是非常重要的。不僅在正常停堆或事故停堆后組件內(nèi)部會出現(xiàn)層流及過渡流動的情況，甚至在反應堆正常運行時的某些組件內(nèi)也會出現(xiàn)。本文將利用CFX計算程序對CEFR帶繞絲棒束組件低Re下的水力特性進行分析。

1 幾何模型

快堆帶繞絲棒束組件內(nèi)的燃料棒通常按照三角形柵格的形式排列，燃料棒之間利用金屬繞絲定位。繞絲按照一定的螺距螺旋式地順時針纏繞在燃料棒上。繞絲的存在大幅增加了幾何建模的難度。在實際結構中，繞絲和燃料棒的接觸為相切接觸，在幾何建模時需將這個點接觸變成面接觸(即將繞絲與燃料棒的接觸由相切變?yōu)橄嘟?，以便進行計算流體域模型建立及網(wǎng)格劃分。

CEFR燃料組件由61根帶繞絲的燃料棒組成，其中組件棒束段的長度為1 350 mm，對整盒組件進行完全模擬需要巨大的計算資源。在有限條件的情況下，本文選用較少螺距的7、19以及61根帶繞絲的燃料棒組成的組件作為對象來進行模擬計算。圖1為1個螺距長度的61根帶繞絲燃料棒的組件，其燃料棒和繞絲的幾何尺寸與CEFR燃料組件的一致：燃料棒直徑Dr=6 mm，繞絲直徑Dw=0.95 mm，燃料棒中心距Pr=7 mm，螺距H=100 mm。

圖1 計算模型示意圖Fig.1 Diagram of calculation model

2 網(wǎng)格劃分

考慮到快堆帶繞絲棒束組件結構的復雜性，不可能采用結構化網(wǎng)格進行網(wǎng)格劃分，本文采用幾何適應性好的非結構化網(wǎng)格，并對繞絲周圍進行加密。

繞絲的存在使組件結構尺寸的尺度變化大，繞絲與相鄰燃料棒之間的間隙非常小，這樣造成網(wǎng)格劃分時非結構化網(wǎng)格的數(shù)量非常龐大，在有限計算資源的情況下很難對組件整體進行計算分析，僅能對很少螺距的帶繞絲棒束組件進行模擬計算。為能獲得組件充分發(fā)展段的水力特性，在僅能對很少螺距的組件進行模擬計算的情況下如果采用常規(guī)進出口條件進行計算，很可能由于入口段的影響，組件內(nèi)的流動達不到充分發(fā)展。因此需選用平移周期性邊界條件來進行處理，利用前一次迭代的出口的結果作為下一次迭代的入口條件，直至達到穩(wěn)定收斂時進出口的結果一致以模擬組件內(nèi)充分發(fā)展的流動。

3 邊界條件

3.1 1個螺距帶繞絲棒束組件

為得到充分發(fā)展段的帶繞絲棒束組件的水力特性，結合繞絲的周期性，在計算中對進出口平面采用平移周期邊界條件，指定通過的質量流量。由于采用了周期邊界條件，僅需計算1個螺距的帶繞絲棒束組件即可，可在有限的計算資源下計算更多數(shù)量的燃料棒組成的組件。

計算中使用液態(tài)金屬鈉作為工質，CEFR燃料組件入口鈉溫度為360 ℃，出口鈉溫度為530 ℃，故選取組件出入口平均溫度445 ℃為定性溫度確定鈉的物性。

3.2 4個螺距帶繞絲棒束組件

通過對4個螺距帶繞絲棒束組件的水力特性進行分析，可得到其入口、出口段的長度以及充分發(fā)展段水力特性變化規(guī)律。4個螺距的帶繞絲棒束組件需要大量的計算資源，本文僅進行了7根燃料棒組成的組件的水力特性分析。

由于在此部分計算中要對組件出入口情況進行分析，故不采用周期性邊界條件。設置入口質量流量，同時出口為開口邊界，邊界壓力設為0，其他條件與3.1節(jié)中一致。

4 水力特性結果

4.1 1個螺距帶繞絲棒束組件

由于繞絲的交混作用，帶繞絲棒束組件在很低Re時已進入過渡工況。本文分別使用層流模型和SST轉捩模型進行計算并將計算結果與國際上公開發(fā)表的經(jīng)驗Engel關系式[1]、Cheng關系式[2]以及劉一哲針對CEFR組件提出的ICRT關系式[3]的計算結果進行比較，從而得到與實際情況比較接近的計算模型。

Engel關系式如下。

層流：

(1)

湍流：

(2)

過渡流動：

400≤Re≤5 000

(3)

式中：f為摩擦阻力系數(shù)；Ψ為間斷因子，是一連接湍流區(qū)和層流區(qū)的參數(shù)，Ψ=(Re-400)/4 600。

Cheng關系式如下。

層流：

(4)

湍流：

(5)

過渡流動：

ReL

(6)

其中：

式中：CfL、CfT分別為層流和湍流關系式的系數(shù)；ReL為層流臨界雷諾數(shù)；ReT為湍流臨界雷諾數(shù)；Pt為棒束柵格距。

ICRT關系式如下。

層流：

(7)

湍流：

f=g1d1fs1M1X1+g2d2fs2M2X2

Re>5 000

(8)

過渡流動：

400≤Re≤5 000

(9)

式中：φ為與水力當量粗糙度相關的參數(shù)；g1、g2為面積分配因子；fs1、fs2為無繞絲通道的表面摩擦因子；d1、d2為棒束通道平均水力直徑與子通道水力直徑的比值；M1、M2為由繞絲引起的摩擦倍數(shù)；X1、X2為不同通道的割流參數(shù)[3]。對于CEFR，φ=25.8，M1=1.124，湍流區(qū)參數(shù)本文不做詳細介紹。

考慮到各類型組件存在不同的棒束數(shù)量，本文計算了61、19及7根燃料棒的組件。

1) 61根帶繞絲棒束組件的水力特性

圖2 61根帶繞絲棒束組件水力特性計算結果與實驗以及各關系式的比較Fig.2 Comparison of hydraulic characteristic calculation result of 61 wire-wrapped rod bundle assembly with experimental result and different relational expressions

61根帶繞絲棒束組件水力特性計算結果如圖2所示。從圖2可看出，利用SST轉捩模型計算的結果在大多情況下與俄羅斯實驗結果[15]符合較好，尤其是在Re為300～2 000時，計算值與實驗值的偏差較小，隨著Re變小或Re變大時這一偏差逐步增大。層流模型計算結果與實驗數(shù)據(jù)以及經(jīng)驗關系式的偏差均較大；而在上述各經(jīng)驗關系式中，Engel關系式的計算結果與實驗結果符合最好，在絕大多數(shù)情況下與實驗值的偏差小于5%；而ICRT關系式的計算結果在Re較小時與實驗值符合較好，而在Re較大時與實驗值的偏差隨著Re的增大而變大；Cheng關系式的計算結果整體與實驗值偏差較大，尤其是在過渡工況的起始Re附近偏差最大。因此，對于CEFR燃料組件這種幾何參數(shù)下的帶繞絲棒束組件，利用Engel關系式來進行低Re下組件的阻力系數(shù)計算準確程度較高。

圖3示出Re=534時與入口距離為0、1/3螺距、2/3螺距、1個螺距平面位置的速度場分布。可看出，由于采用平移周期性邊界條件，出入口的流場保持一致，可模擬充分發(fā)展段帶繞絲棒束組件的流動。同時，1/3螺距、2/3螺距平面位置的速度場分布可通過將入口速度場分別順時針旋轉120°、240°得到，這主要是繞絲順時針旋轉纏繞在燃料棒產(chǎn)生的效果。值得注意的是，由于燃料棒與燃料盒之間的間距較大使得流動阻力較小，造成這些區(qū)域的冷卻劑流速最大，位于中心的燃料棒周圍的流速反而不大，這對于中心位置的燃料棒的冷卻是不利的。

2) 19根帶繞絲棒束組件的水力特性

19根帶繞絲棒束組件計算結果如圖4、5所示。與61根帶繞絲棒束組件的結果相同，利用SST轉捩模型計算得到的結果大多情況下與Engel關系式符合較好(差別小于5%)，與ICRT關系式結果的差別稍大。而層流模型計算結果與經(jīng)驗關系式差別較大。因此利用SST轉捩模型分析繞絲幾何參數(shù)一致的19根燃料棒組件的阻力系數(shù)也將與實際情況符合較好。

3) 7根帶繞絲棒束組件的水力特性

7根帶繞絲棒束組件計算結果如圖6、7所示。將層流模型和SST轉捩模型計算的結果與Engel關系式及ICRT模型結果比較發(fā)現(xiàn)，SST轉捩模型計算的結果在大多情況下與Engel關系式符合較好(差別小于5%)，僅在Re<300之后差別逐漸增大。而該計算結果與ICRT模型結果僅在Re=400附近較接近，隨著Re的增大或減小，這種差別逐漸增大。而利用層流模型計算的結果則與上述兩個關系式相差較大：層流模型計算結果與上述Engel關系式結果的差別保持在36%～48%，與ICRT關系式結果的差別保持在30%～36%。因此對于7根帶繞絲棒束組件的阻力系數(shù)，SST轉捩模型的計算結果也與實際情況符合較好。

圖3 61根帶繞絲棒束組件不同截面速度分布Fig.3 Velocity distribution for different cross sections of 61 wire-wrapped rod bundle assembly

圖4 19根帶繞絲棒束組件水力特性計算結果與各關系式比較Fig.4 Comparison of hydraulic characteristic calculation result of 19 wire-wrapped rod bundle assembly with different relational expressions

綜上所述，通過將61、19及7根燃料棒的帶繞絲棒束組件的水力特性計算結果與實驗以及經(jīng)驗關系式對比發(fā)現(xiàn)，利用SST轉捩模型計算低Re下組件水力特性的結果與實驗及經(jīng)驗關系式符合較好，因此該模型可應用于低Re下此類繞絲結構的分析計算。在后續(xù)的分析計算中均使用該模型。

4.2 4個螺距帶繞絲棒束組件

帶繞絲棒束組件的水力特性對于組件設計、堆芯流量分配以及事故工況下的流動分析都具有重要意義。由于繞絲帶來的橫向流動，帶繞絲棒束組件的水力特性實驗測量時需考慮橫向流動帶來的影響。

由于4個螺距的帶繞絲棒束組件水力特性分析計算需要相當龐大的網(wǎng)格數(shù)量，因此本文選用7根燃料棒的帶繞絲棒束組件進行組件出、入口段及充分發(fā)展段水力特性的分析。

在4個螺距的7根帶繞絲棒束組件水力特性的計算中，計算模型采用SST轉捩模型。將組件盒壁面的6個面編號，分別得到每個面與距離組件入口0、1/4螺距、1/3螺距、1/2螺距、2/3螺距、1個螺距、2個螺距、3個螺距、10/3螺距、11/3螺距、4個螺距的平面的交線中點的靜壓。組件各面編號及長度方向位置(軸向各位置分別與圖中a～k字母對應)如圖8所示。

圖5 19根帶繞絲棒束組件不同截面速度分布Fig.5 Velocity distribution for different cross sections of 19 wire-wrapped rod bundle assembly

圖6 7根帶繞絲棒束組件水力特性計算結果與各關系式比較Fig.6 Comparison of hydraulic characteristic calculation result of 7 wire-wrapped rod bundle assembly with different relational expressions

1) 入口段壓力分布

比較組件6個組件盒壁面在a、c、e、f處(分別對應0、1/3螺距、2/3螺距、1個螺距的位置)的靜壓，如圖9所示。從圖中可看出，在同一軸向位置，組件6個面上的壓力不一致，這也是由于繞絲帶來的橫向流動產(chǎn)生的影響。這種各面壓力不一致的差別會隨著Re的增大而越發(fā)明顯[4]。位于a和f處各面的壓力分布趨勢基本一致，而兩個中間位置c和e的壓力分布有所不同。

實際上，由于繞絲周期性的纏繞，各面的壓力分布也存在著周期性的旋轉變化。對圖9中

圖7 7根帶繞絲棒束組件不同截面速度分布Fig.7 Velocity distribution for different cross sections of 7 wire-wrapped rod bundle assembly

數(shù)據(jù)稍加處理，即將c、e面上6個點與a或f面上6個點根據(jù)繞絲的旋轉方向進行對應(類似“移動相位”)，即c面上點3對應a面上的點1，e面上點5對應a面上點1，以此類推。同時將每個面上各點的壓力換算成與該面各點壓力平均值的差(p-p平面平均)，處理后的結果如圖10所示。從圖10可看出，經(jīng)調整后各壁面的壓力的變化趨勢一致，即由于周期性橫向流動的影響，不同軸向位置的壓力分布與速度分布同樣是旋轉對應的。注意到e和f平面的壓力分布曲線進行相位移動后基本重合，表明e、f位置處的流動已進入充分發(fā)展階段。

圖8 組件壁面和軸向位置Fig.8 Position on wall of assembly and vertical direction

圖9 入口段組件各壁面壓力分布Fig.9 Pressure distribution on different walls at inlet of assembly

2) 出口段壓力分布

出口段壓力分布與入口段的類似，同樣列出各點的壓力分布和進行調整處理后的分布如圖11、12所示。從圖12可看出，經(jīng)調整后各壁面的壓力的變化趨勢一致，即由于周期性橫向流動的影響，不同軸向位置的壓力分布與速度分布同樣是旋轉對應的。注意到各平面的壓力分布曲線進行相位移動后基本重合，表明出口段的影響很小。

圖10 調整后的入口段組件各壁面壓力分布Fig.10 Pressure distribution on different wallsat inlet of assembly after phase change

圖11 出口段組件各壁面壓力分布Fig.11 Pressure distribution on different walls at outlet of assembly

圖12 調整后的出口段組件各壁面壓力分布Fig.12 Pressure distribution on different wallsat outlet of assembly after phase change

3) 整數(shù)倍螺距處壓力分布

圖13為組件上各整數(shù)倍螺距處的壓力分布。從圖13可發(fā)現(xiàn)，從第1個螺距位置的壓力分布開始到第4個螺距位置的壓力分布曲線相互平行，意味著流動進入到充分發(fā)展。相互平行的曲線意味著同一面上1個螺距的壓降均相等，而且與其他面上1個的螺距的壓降一致，并且同一面上整數(shù)倍螺距的壓降恰好是1個螺距壓降的相應整數(shù)倍，這說明雖然繞絲產(chǎn)生的橫向流動使組件6個壁面上壓力分布有所不同，但同一壁面上壓降沿著軸向按螺距均勻分布。這樣在充分發(fā)展段，組件的壓降沿著軸向也是按螺距均勻分布的，組件每個螺距的阻力系數(shù)不變。這一結論也與文獻[16]中的實驗結果相符。在進行帶繞絲棒束組件水力特性測量時，需在組件同一面上按照整數(shù)倍螺距來布置測點，才能避免由于橫向流動對測量帶來的影響。

圖13 整數(shù)倍螺距處組件各壁面壓力分布Fig.13 Pressure distribution on different walls at positions of integer multiples of helix

通過上述結果可知，在較低Re下該帶繞絲棒束組件的入口穩(wěn)定段長度小于1/2螺距的，出口段的影響很小，其余部分為充分發(fā)展段。在充分發(fā)展段內(nèi)部，流體流動按照螺距呈周期性分布，因此，在分析計算該類帶繞絲棒束組件時，只需保證入口、出口段以及至少1個螺距的充分發(fā)展段即可充分說明組件內(nèi)流體的流動情況，其余為周期性部分。因此本文中4個螺距的計算結果能涵蓋全尺寸組件。

5 結論

1) 通過將CFX計算結果與俄羅斯實驗數(shù)據(jù)以及國際上的相關經(jīng)驗關系式比較可得到，利用SST轉捩模型計算帶繞絲棒束組件低Re下的水力特性可得到比較好的結果，因此利用該模型計算低Re下CEFR燃料組件的水力特性是合適的。

2) 通過對4個螺距的7根帶繞絲棒束組件的水力特性分析表明，在較低Re下該帶繞絲棒束組件的入口穩(wěn)定段長度小于1/2螺距的，出口段的影響很小。同時，在流動達到充分發(fā)展后，雖然繞絲產(chǎn)生的橫向流動使組件6個壁面上壓力分布有所不同，但同一壁面上壓降沿著軸向按螺距是均勻分布的，從而組件的壓降沿著軸向也是按螺距均勻分布的，組件每個螺距的阻力系數(shù)不變。

3) 計算結果給實驗的指導：在進行帶繞絲棒束組件水力特性測量時，需在組件同一面上按照整數(shù)倍螺距來布置測點，才能避免由于橫向流動對測量帶來的影響。