秋穗正，張澤秦，張智鵬，王成龍，郭凱倫，田文喜，蘇光輝

(西安交通大學(xué) 核科學(xué)與技術(shù)學(xué)院，陜西 西安 710049)

無人水下潛航器(UUV)是探索海洋奧秘、發(fā)掘海洋資源、維護(hù)海洋安全的有效任務(wù)載體，可執(zhí)行長期無人自主化海洋任務(wù)。UUV的動力源是決定其可執(zhí)行任務(wù)范圍的關(guān)鍵技術(shù)，而傳統(tǒng)水下電源都存在諸如燃料電池質(zhì)量重、體積大、價(jià)格昂貴，柴油機(jī)結(jié)構(gòu)復(fù)雜、噪聲大，蓄電池續(xù)航能力弱、污染重等問題，無法滿足UUV對高能量密度、長續(xù)航、高可靠性電源的需求[1]。與傳統(tǒng)電源相比，核動力電源具有更高的能量密度、更長的使用壽命、更高的可靠性等特點(diǎn)，是UUV電源非常理想的方案。

熱管冷卻核反應(yīng)堆(簡稱熱管反應(yīng)堆)基于固態(tài)反應(yīng)堆設(shè)計(jì)理念，采用液態(tài)金屬高溫?zé)峁芤苑悄軇臃绞綄?dǎo)出堆芯熱量。整個堆芯無需冷卻劑回路，可大幅簡化堆芯設(shè)計(jì)，縮小反應(yīng)堆體積，并提高反應(yīng)堆固有安全性。熱管反應(yīng)堆具有結(jié)構(gòu)緊湊、尺寸小、固有安全性高、模塊化建造和運(yùn)行的技術(shù)特點(diǎn)，能很好地滿足空[2-4]/陸/海[5]等特定應(yīng)用場景、分布式電力供應(yīng)、偏遠(yuǎn)地區(qū)能源供應(yīng)[6]等領(lǐng)域的需求，因此熱管反應(yīng)堆得到了快速發(fā)展，研究者針對熱管反應(yīng)堆的關(guān)鍵技術(shù)開展了廣泛研究，如熱管運(yùn)行特性[7-9]、堆芯多物理場耦合特性[10]、熱管反應(yīng)堆系統(tǒng)運(yùn)行安全特性[11]、熱管與熱電轉(zhuǎn)換裝置耦合特性[12]、熱管反應(yīng)堆原型樣機(jī)制造[13-14]等。

西安交通大學(xué)熱工水力研究室在高溫?zé)峁芗盁峁芊磻?yīng)堆領(lǐng)域耕耘近10年，在高溫?zé)峁軅鳠釋?shí)驗(yàn)[15-16]、高溫?zé)峁苓\(yùn)行數(shù)值模擬[17-19]、大功率高溫?zé)峁茉O(shè)計(jì)制造[20]、固態(tài)熱管反應(yīng)堆設(shè)計(jì)[21]、熱管反應(yīng)堆系統(tǒng)安全分析[22]等研究方向開展了實(shí)驗(yàn)研究和理論分析。自2019年12月起，由西安交通大學(xué)牽頭，依托國家重點(diǎn)研發(fā)計(jì)劃項(xiàng)目“新型海洋靜默式熱管核反應(yīng)堆技術(shù)研究”，集合國內(nèi)優(yōu)勢單位，著力突破海洋核動力UUV“卡脖子”關(guān)鍵技術(shù)，為我國熱管反應(yīng)堆工程化應(yīng)用奠定設(shè)計(jì)基礎(chǔ)和理論依據(jù)。本文面向海洋重型UUV的動力需求，提出新型海洋靜默式熱管反應(yīng)堆NUSTER-100設(shè)計(jì)方案，建立反應(yīng)堆全系統(tǒng)數(shù)學(xué)物理模型，并開發(fā)先進(jìn)熱管反應(yīng)堆瞬態(tài)分析程序HEART，深入研究NUSTER-100穩(wěn)態(tài)、冷啟動瞬態(tài)及反應(yīng)性引入瞬態(tài)工況下的熱工水力及安全特性，為熱管應(yīng)用于熱管反應(yīng)堆及其構(gòu)成的核能電源系統(tǒng)提供設(shè)計(jì)及理論支撐。

1 海洋靜默式熱管反應(yīng)堆設(shè)計(jì)

1.1 反應(yīng)堆整體結(jié)構(gòu)

NUSTER-100整體結(jié)構(gòu)如圖1所示，主要包括堆芯及屏蔽系統(tǒng)、能量轉(zhuǎn)換系統(tǒng)、余熱排出系統(tǒng)及堆芯應(yīng)急冷卻系統(tǒng)。堆芯產(chǎn)生的熱量由熱管傳輸至溫差發(fā)電模塊，發(fā)電后剩余的廢熱由冷卻水板及余熱排出系統(tǒng)傳輸至冷卻水箱[23]。整個系統(tǒng)無需轉(zhuǎn)動部件，反應(yīng)堆一回路沒有冷卻劑流動。NUSTER-100系統(tǒng)布置方式如圖2所示，系統(tǒng)主要設(shè)計(jì)參數(shù)列于表1。

表1 NUSTER-100系統(tǒng)主要設(shè)計(jì)參數(shù)

圖1 NUSTER-100系統(tǒng)結(jié)構(gòu)示意圖

圖2 NUSTER-100系統(tǒng)布置示意圖

1.2 堆芯布置

NUSTER-100堆芯采用11×11單元式布置，共109個堆芯單元。每個單元以熱管為中心，周圍布置12根燃料棒。堆芯共有燃料棒480根，采取72%/50%/19.75%富集度UO2三區(qū)布置方式；燃料棒與熱管之間采用Mo作為填充基體，整個反應(yīng)堆外層包圍有3.3 mm厚的絕熱層。反應(yīng)堆堆芯和熱管單元如圖3所示，堆芯及熱管單元主要設(shè)計(jì)參數(shù)列于表2。

圖3 反應(yīng)堆堆芯和熱管單元示意圖

表2 堆芯及熱管單元主要設(shè)計(jì)參數(shù)

1.3 能量轉(zhuǎn)換系統(tǒng)

熱管反應(yīng)堆采用高溫鈉熱管作為熱量傳輸裝置，采用靜默式熱電轉(zhuǎn)換作為能量轉(zhuǎn)換方式。為預(yù)留空間給反應(yīng)性控制部件的驅(qū)動機(jī)構(gòu)，堆芯中109根熱管呈非對稱布置，其中一側(cè)布置58根熱管，共6層；另一側(cè)布置51根，共5層。每層熱管由高熱導(dǎo)率材料制造的固體換熱器箍緊限位，固體換熱器上布置一層熱電轉(zhuǎn)換單元串聯(lián)構(gòu)成的發(fā)電矩陣。熱電轉(zhuǎn)換單元的冷端為冷卻水板，為整個能量轉(zhuǎn)換系統(tǒng)提供冷源。熱管-熱電轉(zhuǎn)換器件-冷卻水板構(gòu)成“三明治”結(jié)構(gòu)，兩側(cè)內(nèi)部的冷卻水板被兩個“三明治”結(jié)構(gòu)共用。能量轉(zhuǎn)換系統(tǒng)設(shè)計(jì)如圖4a所示，熱電轉(zhuǎn)換單元結(jié)構(gòu)如圖4b所示，熱電器件的主要設(shè)計(jì)參數(shù)列于表3。

表3 熱電轉(zhuǎn)換單元的主要設(shè)計(jì)參數(shù)

a——能量轉(zhuǎn)換系統(tǒng)；b——熱電轉(zhuǎn)換單元

2 熱管反應(yīng)堆數(shù)學(xué)物理模型

針對海洋靜默式熱管反應(yīng)堆設(shè)計(jì)特點(diǎn)，建立包括堆芯功率模型、堆芯通道傳熱模型、熱管傳熱模型、熱電轉(zhuǎn)換模型及冷端換熱模型等在內(nèi)的熱管反應(yīng)堆系統(tǒng)分析模型，實(shí)現(xiàn)對NUSTER-100全系統(tǒng)穩(wěn)態(tài)和瞬態(tài)熱工水力特性的分析。

2.1 堆芯功率模型

采用考慮多組緩發(fā)中子的點(diǎn)堆中子動力學(xué)方程求解堆芯裂變功率，并考慮各類反饋效應(yīng)引入的反應(yīng)性。

點(diǎn)堆動力學(xué)方程表示如下：

(1)

(2)

式中：n為堆內(nèi)中子密度，s-1·m-3；t為時(shí)間，s；ρ為總反應(yīng)性；β為總的有效緩發(fā)中子份額；Λ為中子代時(shí)間，s；λi為第i組緩發(fā)中子的衰變常量，s-1；Ci為第i組緩發(fā)中子先驅(qū)核的濃度，m-3；βi為第i組緩發(fā)中子份額；nc為緩發(fā)中子組數(shù)。

點(diǎn)堆方程中，堆芯裂變功率主要由總反應(yīng)性ρ控制，任一時(shí)刻總的反應(yīng)性可表示為：

ρ(t)=ρCR(t)+ρDOP(t)+ρEXPAN(t)+ρHP(t)

(3)

式中：ρ(t)為總反應(yīng)性；ρCR(t)為控制棒和停堆棒等引入的反應(yīng)性；ρDOP(t)為燃料的多普勒反饋；ρEXPAN(t)為反射層與基體膨脹引入的負(fù)反饋；ρHP(t)為熱管元件材料與其內(nèi)部工質(zhì)因溫度與物性變化引入的反應(yīng)性反饋。

2.2 堆芯通道傳熱模型

熱管反應(yīng)堆的堆芯單元結(jié)構(gòu)復(fù)雜，從燃料芯塊到熱管蒸氣空間包含6～7層固體材料區(qū)域。對于單個堆芯單元，由于燃料棒與熱管的布置并非完全對稱，這些固體區(qū)域的熱特性存在差異；但對于堆芯整體，這種差異是可忽略的。因此，為簡化建模，可等效地將熱管單元內(nèi)不規(guī)則形狀的區(qū)域視為規(guī)則的圓環(huán)，矩形堆芯通道可按照近似圓柱形通道進(jìn)行建模[22]，堆芯單元的等效示意圖如圖5所示。每根燃料棒按照真實(shí)的圓柱導(dǎo)熱模型計(jì)算，包殼、基體區(qū)域按照等效圓環(huán)模型計(jì)算，熱管區(qū)域則與熱管計(jì)算模塊進(jìn)行邊界條件傳遞。

圖5 堆芯單元等效示意圖

根據(jù)等效的圓柱形堆芯單元以及燃料芯塊導(dǎo)熱的特點(diǎn)，可忽略各層結(jié)構(gòu)材料中的內(nèi)熱源，認(rèn)為內(nèi)熱源只存在于燃料芯塊中，則芯塊區(qū)域的傳熱為有內(nèi)熱源的導(dǎo)熱，控制方程如下：

(4)

式中：ρU為燃料芯塊的密度，kg/m3；cU為燃料芯塊的比熱，J/(kg·K)；TU為燃料芯塊的溫度，K；λU為燃料芯塊的導(dǎo)熱系數(shù)，W/(m·K)；QV為燃料控制體的熱源密度，W/m3。

由于燃料包殼與燃料的基體材料均為Mo，可視為同一種導(dǎo)熱基體，其傳熱均為無內(nèi)熱源的導(dǎo)熱，采用統(tǒng)一的控制方程：

(5)

式中：ρM為導(dǎo)熱基體的密度，kg/m3；cM為導(dǎo)熱基體的比熱，J/(kg·K)；TM為導(dǎo)熱基體的溫度，K；λM為導(dǎo)熱基體的導(dǎo)熱系數(shù)，W/(m·K)。

由于熱導(dǎo)問題的特殊性，上述圓環(huán)等效建模一般有兩種方式：第1種是保證導(dǎo)熱系數(shù)計(jì)算的真實(shí)性，即各等效圓環(huán)的厚度采用真實(shí)傳熱特征長度，為達(dá)到計(jì)算的準(zhǔn)確性，需要對熱容進(jìn)行等效計(jì)算；第2種是保證熱容計(jì)算的真實(shí)性，即認(rèn)為等效圓環(huán)的體積與原區(qū)域體積相同，則傳熱長度增大，需要對導(dǎo)熱系數(shù)進(jìn)行等效計(jì)算。第1種等效方式在處理穩(wěn)態(tài)計(jì)算問題上十分奏效，但對于瞬態(tài)計(jì)算偏差較大；第2種等效方式適用性更好，對導(dǎo)熱系數(shù)的處理也較簡單方便，其計(jì)算正確性已被證明[24]。對于包殼和基體區(qū)域的等效模型，采用第2種建模方案，則各區(qū)域的等效厚度可采用下式計(jì)算：

(6)

式中：δn,e為第n層等效圓環(huán)的厚度，m；An為第n層等效圓環(huán)的真實(shí)面積，m2；rn-1,o為第n-1層等效圓環(huán)的半徑，m；下標(biāo)e代表等效值，o代表原始值。

各層等效圓環(huán)的修正導(dǎo)熱系數(shù)按下式計(jì)算[25]：

(7)

式中：λn,e和λn,r分別為第n層等效圓環(huán)的等效導(dǎo)熱系數(shù)和真實(shí)導(dǎo)熱系數(shù)，W/(m·K)；rn,e為第n層等效圓環(huán)的半徑，m；rn-1,e為第n-1層等效圓環(huán)的半徑，m；rn,r為第n層等效圓環(huán)的真實(shí)傳熱特征長度，m，對于基體取加權(quán)平均值；下標(biāo)r代表真實(shí)值。

燃料區(qū)域與包殼間的氦氣氣隙在計(jì)算中視為邊界條件進(jìn)行處理。氦氣氣隙的熱容非常小，因此可直接采用上述第1種計(jì)算方式，即認(rèn)為計(jì)算中的氣隙寬度等于氦氣氣隙的真實(shí)傳熱特征長度。則燃料計(jì)算域與等效計(jì)算域間的邊界條件為：

(8)

式中：λHe為氦氣的導(dǎo)熱系數(shù)，W/(m·K)；δHe為氦氣氣隙寬度，m；εUM為燃料芯塊與包殼材料間的表面發(fā)射率；σBoz為斯蒂芬-玻耳茲曼常數(shù)。

2.3 熱管傳熱模型

熱管反應(yīng)堆采用堿金屬高溫?zé)峁軅鬟f堆芯熱量，其傳熱模型如圖6所示。在軸向上，熱管一般劃分為蒸發(fā)段、絕熱段和冷凝段，分別影響模型的外邊界條件；在徑向上，熱管一般劃分為管壁、吸液芯和蒸氣空間，分別影響模型的控制方程。

圖6 熱管的傳熱模型示意圖

高溫?zé)峁茉诠ぷ鲿r(shí)，吸液芯內(nèi)的堿金屬工質(zhì)因熱導(dǎo)率較高而流速較低，可忽略吸液芯內(nèi)液態(tài)工質(zhì)的流動過程，將吸液芯區(qū)域和管壁區(qū)域都按照純導(dǎo)熱進(jìn)行計(jì)算；另一方面，熱管要將熱量導(dǎo)出堆芯，其沿軸向的熱量傳遞不能忽略，因此熱管的管壁和吸液芯區(qū)域采用二維導(dǎo)熱控制方程進(jìn)行建模。

熱管吸液芯區(qū)域被視為液態(tài)工質(zhì)與固態(tài)材料的混合固態(tài)基體，該基體的體積熱容和導(dǎo)熱系數(shù)由其組成部分共同決定。對于蒸氣區(qū)域，根據(jù)蒸氣連續(xù)流態(tài)是否完全建立，采取不同控制方程[19]。不同區(qū)域的計(jì)算模型如圖7所示。

圖7 熱管傳熱計(jì)算模型示意圖

2.4 熱電轉(zhuǎn)換模型

熱管反應(yīng)堆系統(tǒng)一般采用靜態(tài)發(fā)電策略，如熱離子發(fā)電、溫差發(fā)電等，也有部分采用動態(tài)發(fā)電策略，如斯特林發(fā)電等。本文提出的NUSTER-100采用靜默式溫差發(fā)電實(shí)現(xiàn)熱電轉(zhuǎn)換。系統(tǒng)中每個熱管單元兩端各對應(yīng)1組溫差發(fā)電器件。在熱管冷凝段，不同軸向位置上的熱電轉(zhuǎn)換計(jì)算簡化為沿?zé)峁軓较蚍较虻囊痪S計(jì)算，其熱端邊界為熱管冷凝段外表面，冷端邊界為冷卻水板外表面，中間依次經(jīng)歷集熱板、銅極墊片、熱電轉(zhuǎn)換單元、銅極墊片，如圖8所示。由于熱管冷凝段外圍與高熱導(dǎo)率集熱板相連接，可認(rèn)為溫差發(fā)電器件的熱端溫度相同。

圖8 熱電轉(zhuǎn)換單元示意圖

上述單元中，集熱板和銅極墊片中發(fā)生的都是純導(dǎo)熱過程，其導(dǎo)熱方程分別為：

(9)

(10)

式中：ρHCP為集熱板材料的密度，kg/m3；cHCP為集熱板材料的比熱，J/(kg·K)；THCP為集熱板材料的溫度，K；λHCP為集熱板材料的導(dǎo)熱系數(shù)，W/(m·K)；ρCu為銅的密度，kg/m3；cCu為銅的比熱，J/(kg·K)；TCu為銅的溫度，K；λCu為銅的導(dǎo)熱系數(shù)，W/(m·K)。

熱電轉(zhuǎn)換單元的節(jié)點(diǎn)溫度方程為：

(1-η)-ACqC(TTEG-TC)

(11)

式中：η為熱電轉(zhuǎn)換效率；ρTEG為熱電轉(zhuǎn)換單元的密度，kg/m3；cTEG為熱電轉(zhuǎn)換單元的比熱，J/(kg·K)；q為熱流密度，W/m2；A為垂直于熱流方向的橫截面積，m2；下標(biāo)H和C分別代表熱電轉(zhuǎn)換單元的熱端和冷端。

熱電轉(zhuǎn)換效率η的求解綜合考慮了塞貝克效應(yīng)、珀?duì)柼?yīng)和湯姆遜效應(yīng)[26]，由下式計(jì)算：

η=

(12)

電流I由下式求解：

(13)

(14)

式中：V為電壓；αeff為等效塞貝克系數(shù)，V/K；αp,i為P型電偶臂第i段材料的塞貝克系數(shù)，V/K；αn,i為N型電偶臂第i段材料的塞貝克系數(shù)，V/K；ΔTp,i為P型電偶臂第i段材料的溫差，K；ΔTn,i為N型電偶臂第i段材料的溫差，K；m為熱電轉(zhuǎn)換單元總段數(shù)；RL為外接負(fù)載值，Ω；λeff,p和λeff,n分別為熱電轉(zhuǎn)換單元P型和N型電偶臂的等效導(dǎo)熱系數(shù)；Rp和Rn分別為P型和N型電偶臂的等效電阻。

(15)

(16)

式中：λ為溫差發(fā)電材料的導(dǎo)熱系數(shù)，W/(m·K)；A為垂直于熱流方向的橫截面積，m2；L為單腿長度，m；σ為溫差發(fā)電材料電導(dǎo)，S/m；下標(biāo)x為n時(shí)代表N型電偶臂，為p時(shí)代表P型電偶臂，i代表熱電材料層數(shù)。

上述求解過程中，冷端溫度TC與熱電轉(zhuǎn)換效率η存在耦合關(guān)系，因此需要迭代求解。由于負(fù)載RL確定后，η僅是溫差(TH-TC)的函數(shù)，首先建立溫差(TH-TC)與η的關(guān)系η=f(TH-TC)，假定熱電轉(zhuǎn)換效率為0，求出冷端溫度TC，再由函數(shù)關(guān)系f得到迭代初值。為節(jié)省計(jì)算資源，可直接建立η與TC和TH-TC的二元關(guān)系η=F[(TH-TC)，TC]，通過插值進(jìn)行求解。

2.5 冷端換熱模型

熱管反應(yīng)堆系統(tǒng)中，冷卻水板與熱電轉(zhuǎn)換元件的冷端進(jìn)行換熱，將發(fā)電后剩余熱量排出堆芯。冷卻水板內(nèi)部設(shè)有矩形流道，由外部系統(tǒng)提供恒定入口壓力、入口溫度的冷卻水，系統(tǒng)最終冷阱為海水。冷卻水板內(nèi)部的冷卻水建模為一維不可壓縮流動，其質(zhì)量、動量和能量守恒方程為：

(17)

(18)

(19)

式中：ρ為冷卻劑密度，kg/m3；W為冷卻劑質(zhì)量流量，kg/s；A為冷卻水流道截面積，m2；p為冷卻劑壓力，Pa；h為冷卻劑焓值，kJ/kg；De為冷卻水流道當(dāng)量直徑，m；qw為傳入冷卻水的總熱流密度，W/m2；Π為冷卻水流道的濕周，m；S為沿流動方向的軸向位置，m。

3 數(shù)值算法與程序開發(fā)驗(yàn)證

3.1 數(shù)值算法

通過建立熱管反應(yīng)堆全系統(tǒng)數(shù)學(xué)物理模型，構(gòu)建了相應(yīng)的求解熱管反應(yīng)堆系統(tǒng)熱工水力特性的封閉控制方程組，可采用多種數(shù)值算法進(jìn)行求解。

上述熱管反應(yīng)堆系統(tǒng)控制方程組的基本形式為熱工水力參數(shù)對時(shí)間的導(dǎo)數(shù)與參數(shù)本身之間的顯式或隱式函數(shù)關(guān)系，通過離散可轉(zhuǎn)化為非線性常微分方程組的初值問題，其具有如下形式：

(20)

式中：y0為給定的初值；f為由控制方程得到的函數(shù)關(guān)系。

Gear算法[27]采用隱式向后差分方法，設(shè)計(jì)了一種針對剛性方程組的穩(wěn)定策略。Gear算法采用牛頓迭代法進(jìn)行隱式求解，利用矩陣的系數(shù)結(jié)構(gòu)特點(diǎn)直接求解線性方程，因此每前進(jìn)一個步長解隱式方程組所需要的工作量較小，從而加快了計(jì)算速度。本研究采用Gear算法求解非線性常微分方程組的初值問題。

熱管反應(yīng)堆系統(tǒng)控制方程組中除非線性常微分初值問題外，還存在非線性常微分邊值問題，其與方程(20)具有相同的構(gòu)造形式，如熱管蒸氣區(qū)內(nèi)速度、壓力、密度等參數(shù)構(gòu)成的一階非線性方程組等。針對這類問題，采用顯式Runge-Kutta法(R-K法)進(jìn)行求解。

3.2 程序開發(fā)架構(gòu)

在上述熱管反應(yīng)堆系統(tǒng)數(shù)學(xué)物理模型與數(shù)值算法的基礎(chǔ)上，開發(fā)了具有自主知識產(chǎn)權(quán)的先進(jìn)熱管反應(yīng)堆瞬態(tài)分析程序HEART。程序采用面向?qū)ο蟮木幊谭妒?，基于模塊化的編程思想，按類化進(jìn)行編制，具有高度的通用性、可讀性、可拓展性與可維護(hù)性。各系統(tǒng)均可調(diào)用輸入模塊、初始化模塊和控制方程計(jì)算模塊，各模塊既可獨(dú)立運(yùn)行，又可作為子模塊由主程序調(diào)用共同求解。HEART程序總體模塊化設(shè)計(jì)及計(jì)算流程如圖9所示，主要分為物理計(jì)算模塊、熱工計(jì)算模塊以及公用模塊。每個部件模塊內(nèi)部包含每個部件的輸入、初始化、導(dǎo)數(shù)計(jì)算以及邊界傳入傳出部分，刪除或添加一個模塊不會影響其他模塊，各模塊可并行計(jì)算，提升了程序的計(jì)算效率。

圖9 HEART程序結(jié)構(gòu)

3.3 程序模型驗(yàn)證

程序模型的驗(yàn)證與確認(rèn)(V&V)是程序開發(fā)中的重要一環(huán)，對提高程序計(jì)算的可靠性和準(zhǔn)確性具有重要意義。針對熱管反應(yīng)堆，目前國際上公開發(fā)表或內(nèi)部的實(shí)驗(yàn)數(shù)據(jù)較少，因此HEART程序的模型驗(yàn)證與確認(rèn)采用單模塊驗(yàn)證方法。堆芯通道傳熱模型采用商業(yè)程序計(jì)算得到的設(shè)計(jì)值進(jìn)行驗(yàn)證，熱管和熱電轉(zhuǎn)換模型則采用實(shí)驗(yàn)值進(jìn)行驗(yàn)證。

1)堆芯通道傳熱模型驗(yàn)證

使用HEART程序?qū)USTER-100進(jìn)行模擬計(jì)算，軸向功率分布因子由MCNP計(jì)算，結(jié)果如圖10所示。堆芯最熱通道沿軸向和徑向的HEART程序溫度計(jì)算結(jié)果與設(shè)計(jì)值的對比如圖11所示?？煽闯?，程序計(jì)算值與設(shè)計(jì)值的最大偏差不超過5 K，軸向溫度分布計(jì)算值與設(shè)計(jì)值基本一致，證明了堆芯通道傳熱模型建立的合理性與準(zhǔn)確性。

圖10 堆芯軸向功率因子分布

a——沿燃料棒區(qū)域的軸向溫度分布；b——軸向溫度最高點(diǎn)處的徑向溫度分布

2)熱管傳熱模型驗(yàn)證

熱管傳熱模型通過與熱管冷態(tài)啟動實(shí)驗(yàn)進(jìn)行對比驗(yàn)證。實(shí)驗(yàn)所用的熱管詳細(xì)參數(shù)和邊界條件列于表4。HEART程序計(jì)算得到的熱管外壁面沿軸向的溫度分布與實(shí)驗(yàn)值的比較如圖12所示。圖12表明，HEART程序?qū)峁芩矐B(tài)過程的計(jì)算結(jié)果與實(shí)驗(yàn)值吻合較好，在準(zhǔn)穩(wěn)態(tài)部分(25 min后)的計(jì)算偏差小于10 K。

圖12 HEART程序熱管傳熱模型計(jì)算結(jié)果與實(shí)驗(yàn)值的對比

表4 實(shí)驗(yàn)熱管參數(shù)

3)熱電轉(zhuǎn)換模型驗(yàn)證

熱電轉(zhuǎn)換模型通過與單段式溫差發(fā)電器件加熱實(shí)驗(yàn)進(jìn)行對比驗(yàn)證。實(shí)驗(yàn)所用的溫差發(fā)電器件詳細(xì)參數(shù)列于表5。HEART程序計(jì)算得到的溫差發(fā)電器件冷熱端溫度分布與實(shí)驗(yàn)值的對比如圖13所示?？煽闯?，HEART計(jì)算值與實(shí)驗(yàn)所得溫差發(fā)電器件冷熱端溫度變化趨勢吻合較好。計(jì)算值稍高于實(shí)驗(yàn)值的原因是計(jì)算中沒有考慮溫差發(fā)電器件在與冷熱源接觸面上的漏熱和熱阻。功率級別較低的情況下，冷端有足夠的換熱能力，溫度變化很小，僅為2 K左右。計(jì)算結(jié)果的最大偏差小于4 K，證明了熱電轉(zhuǎn)換模型的合理性。

表5 單段式溫差發(fā)電器件參數(shù)

圖13 HEART程序熱電轉(zhuǎn)換模型計(jì)算結(jié)果與實(shí)驗(yàn)值的對比

4 NUSTER-100穩(wěn)態(tài)及瞬態(tài)工況分析

本文在上述熱管反應(yīng)堆數(shù)學(xué)物理模型和HEART程序的基礎(chǔ)上，對NUSTER-100穩(wěn)態(tài)及瞬態(tài)工況下的熱工水力特性進(jìn)行計(jì)算分析。

4.1 滿功率穩(wěn)態(tài)工況分析

NUSTER-100計(jì)算中，HEART程序所需要的點(diǎn)堆中子動力學(xué)參數(shù)和堆芯歸一化功率分布因子均由半確定論蒙特卡羅程序NECP-MCX給出[28]，圖14為NUSTER-100的1/8堆芯通道級絕對功率分布。堆芯外邊界(3.3 mm保溫層)視為絕熱邊界，堆芯的裂變熱全部通過109根熱管傳遞。被多個通道公用的燃料棒由其所在通道的平均燃料棒功率加權(quán)得到，邊燃料棒的功率權(quán)值為1/2，角燃料棒的功率權(quán)值為1/4。

圖14 1/8堆芯通道級絕對功率分布

計(jì)算獲得的NUSTE-100穩(wěn)態(tài)工況下1/8堆芯內(nèi)燃料棒的中心溫度分布如圖15所示。在軸向功率最高點(diǎn)處，中心通道與邊通道的燃料棒中心溫度差為80 K左右，軸向中點(diǎn)處的溫度差為40 K左右。中心通道處燃料棒的功率最高點(diǎn)處與軸向中點(diǎn)處的溫度差約為38 K，與設(shè)計(jì)值一致，而邊通道約為30 K。此結(jié)果表明，熱管反應(yīng)堆的固態(tài)堆芯具有較好的溫度自展平能力，即使堆芯中心通道的功率與外圍通道相差較大，其溫差相對不大，節(jié)省了空間。

a——軸向功率最高點(diǎn)處；b——軸向中點(diǎn)處

1/8堆芯內(nèi)燃料棒的外邊界溫度分布如圖16所示。由于采用了等效模型，針對各通道內(nèi)不同燃料棒的外邊界進(jìn)行了歸一化處理，由此轉(zhuǎn)變?yōu)獒槍νǖ赖臏囟惹蠼?。在燃料棒外邊界處，由于軸向功率分布因子造成的溫度不均勻性已得到緩解。在軸向上，功率最高點(diǎn)和中點(diǎn)處的外邊界溫差下降到11 K左右，這主要得益于熱管和金屬基體的高導(dǎo)熱能力。

a——軸向功率最高點(diǎn)處；b——軸向中點(diǎn)處

1/8堆芯內(nèi)熱管的蒸氣區(qū)溫度分布如圖17所示。由于熱管良好的等溫性，熱管蒸氣區(qū)沿軸向的蒸氣溫降非常低，功率最高點(diǎn)與軸向中點(diǎn)的溫差不到2 K。由于熱管特殊的傳熱性質(zhì)，其蒸發(fā)段溫度與傳遞功率并非呈線性關(guān)系。只要保證堆芯熱管冷凝段的邊界條件一致，堆芯熱管的蒸氣溫度分布便不會發(fā)生明顯變化。熱管的這一特性提供了熱管反應(yīng)堆堆芯溫度自展平能力，并保證了熱電轉(zhuǎn)換模塊工作的穩(wěn)定性。

a——軸向功率最高點(diǎn)處；b——軸向中點(diǎn)處

以中心通道為研究對象，其從燃料棒中心到熱管蒸氣區(qū)的溫度分布如圖18所示。溫度分布不均勻主要發(fā)生在燃料棒區(qū)域，而包殼-基體區(qū)域的溫度分布差異已不明顯，在熱管壁-吸液芯-蒸氣區(qū)域，由于熱管的高導(dǎo)熱能力，沿軸向基本沒有溫差。由于溫度分布的差異較為明顯、不同區(qū)域材料的膨脹系數(shù)各不相同，再加上固體堆芯特殊的裝配方式，熱管反應(yīng)堆堆芯熱-力耦合分析具有必要性。

圖18 中心通道內(nèi)溫度分布

能量轉(zhuǎn)換系統(tǒng)的穩(wěn)態(tài)特性如圖19所示。以堆芯中心通道為研究對象，研究了熱管-固體換熱器-溫差發(fā)電模塊-冷卻劑系統(tǒng)的整體特性，并在熱管蒸發(fā)段采用加密網(wǎng)格，更為準(zhǔn)確地模擬了軸向功率不均勻分布的影響。由圖19可知，蒸發(fā)段與絕熱段間的溫差為40 K，絕熱段與冷凝段間的溫差為38 K。固體換熱器和溫差發(fā)電器件兩端的銅極墊片溫降約為65 K，溫差發(fā)電模塊的平均溫降約為724 K。冷卻水板內(nèi)的冷卻劑吸收余熱后升溫約20 K。假定系統(tǒng)負(fù)載RL與內(nèi)阻(Rn+Rp)相等，則可得到最大電功率。在這種情況下，中心通道的能量轉(zhuǎn)換系統(tǒng)可產(chǎn)生193.5 V開路電壓和1 207.8 W電能，熱電轉(zhuǎn)換效率約為10.88%，具體熱電參數(shù)列于表6。

表6 中心通道能量系統(tǒng)熱電參數(shù)計(jì)算結(jié)果

圖19 能量轉(zhuǎn)換系統(tǒng)溫度分布

4.2 冷態(tài)啟動工況分析

針對熱管反應(yīng)堆啟動方案，目前國內(nèi)外相關(guān)研究較少。多數(shù)研究側(cè)重于穩(wěn)態(tài)分析，或?qū)⒎磻?yīng)堆設(shè)計(jì)為熱態(tài)啟動，即認(rèn)為反應(yīng)堆從一個相對較高的初始溫度進(jìn)行啟動(一般高于所用熱管的工質(zhì)蒸氣轉(zhuǎn)變溫度)，啟堆時(shí)熱管已處于工作狀態(tài)。而對于海洋應(yīng)用環(huán)境，完全冷態(tài)啟動工況的模擬分析是必要的。本文使用HEART程序?qū)USTER-100進(jìn)行冷態(tài)啟動工況分析。

由于NUSTER-100的最熱通道為中心通道，反應(yīng)堆安全啟動的前提是保證中心通道順利啟動。因此選取中心通道作為研究對象，認(rèn)為啟動過程為功率控制，分別選取啟動速率0.1%FP/s、0.05%FP/s和0.02%FP/s進(jìn)行模擬計(jì)算，熱管蒸發(fā)段外壁面最高溫度點(diǎn)如圖20所示?？煽吹?，以0.1%FP/s作為啟動速率時(shí)，在780 s左右熱管遇到毛細(xì)極限，蒸發(fā)段前段溫度飛升230 K左右，熱管發(fā)生啟動失敗。而選取0.05%FP/s以及0.02%FP/s作為啟動速率時(shí)，熱管均可正常啟動。因此，熱管反應(yīng)堆啟動過程中存在啟動速度限值，過快的升功率速度將導(dǎo)致反應(yīng)堆啟動失敗。

圖20 不同啟動速率下中心通道熱管蒸發(fā)段外壁面最高溫度

上述不同啟動速率下中心通道熱管的軸向溫度分布及傳遞功率如圖21所示。由圖21d～f可知，熱管的完全冷態(tài)啟動需要經(jīng)歷功率緩慢增加的過程。在啟動最初的一段時(shí)間內(nèi)，熱管并無功率傳遞，結(jié)合圖21a～c可知，這段時(shí)間內(nèi)熱量主要用于解凍熱管蒸發(fā)段，為熱管的工作死區(qū)。隨著啟動速率的降低，死區(qū)時(shí)間占比逐漸下降，熱管傳遞功率緩慢增加，在這個過程中熱管傳熱主要受聲速極限的制約。聲速極限會導(dǎo)致熱管傳熱值被限制，但一般不會導(dǎo)致熱管失效，結(jié)合圖21a～c可知，這個過程中熱管的主要現(xiàn)象是冷凝段熔化并啟動，連續(xù)流態(tài)蒸氣在整個蒸氣區(qū)建立。突破聲速極限后，熱管功率將迅速增加至與通道功率一致。如果采用線性功率變化，在達(dá)到滿功率水平前，熱管即完成啟動。在熱管能正常啟動的工況下，不論啟動速率如何，熱管啟動所需要的時(shí)間基本為總時(shí)間的75%左右，因此降低熱管啟動速率的本質(zhì)是降低熱管啟動效率來換取安全性。

啟動速率：a，d——0.1%FP/s；b，e——0.05%FP/s；c，f——0.02%FP/s

基于上述研究，為提升熱管反應(yīng)堆的啟動效率，滿足海洋熱管反應(yīng)堆應(yīng)用需求下快速啟動的需求，初步提出了一種三段式熱管反應(yīng)堆啟動方案，如圖22所示。該方案的主要步驟為：1)在啟動初期，采取相對較高的升功率策略，如以0.2%FP/s啟動速率快速解凍熱管蒸發(fā)段，縮短啟動死區(qū)時(shí)間，直到40%FP或熱管升溫速度大幅下降；2)在啟動中期，采取耦合控制升功率策略，監(jiān)測熱管升溫速度，根據(jù)升溫速度采取相對較低的升功率策略，如0.02%FP/s～0.04%FP/s，直到監(jiān)測到冷凝段末端溫度快速升高，認(rèn)為首次突破聲速極限；3)在啟動末期，采取適中的升功率策略，如0.15%FP/s，保證熱管順利進(jìn)入準(zhǔn)穩(wěn)態(tài)階段。

圖22 三段式熱管反應(yīng)堆啟動方案

采用三段式熱管反應(yīng)堆啟動方案重新對NUSTER-100進(jìn)行冷態(tài)啟動計(jì)算，得到的通道溫度變化和熱管傳遞功率如圖23所示。由圖23可知，采用三段式啟動控制方案后，在同樣的啟動時(shí)間內(nèi)，熱管反應(yīng)堆能順利啟動，但啟動過程中熱管不再提前達(dá)到通道功率。在熱管反應(yīng)堆功率處于平臺階段的過程中，堆芯內(nèi)燃料棒中心、包殼表面和熱管表面間的溫差逐漸增加，但增量不大。從總體上看，這種啟動方案有效提高了熱管反應(yīng)堆的啟動效率，保障了啟動過程的安全性。

圖23 三段式熱管反應(yīng)堆啟動方案下的通道溫度和熱管傳遞功率

4.3 反應(yīng)性引入瞬變工況分析

反應(yīng)堆的固有安全性主要來源于堆芯的自穩(wěn)特性。相較于傳統(tǒng)水堆的強(qiáng)負(fù)反饋特性，熱管反應(yīng)堆在宏觀上屬于弱負(fù)反饋型反應(yīng)堆，需要根據(jù)熱管反應(yīng)堆的設(shè)計(jì)方案對其自穩(wěn)能力進(jìn)行評估。本文模擬計(jì)算NUSTER-100的反應(yīng)性引入瞬變工況。

圖24為+10 pcm反應(yīng)性引入工況下反應(yīng)堆核功率、堆芯最高溫度、總反應(yīng)性的瞬態(tài)響應(yīng)過程。由圖24可知，滿功率下，+10 pcm反應(yīng)性引入后，反應(yīng)堆系統(tǒng)需要350 s左右達(dá)到新的穩(wěn)態(tài)，新穩(wěn)態(tài)下堆芯核功率為1 019 kW左右，穩(wěn)態(tài)功率上升了1.9%，核功率的超調(diào)量達(dá)到4.85%，堆芯最高溫度上升近30 K。此外，由于金屬基體的熱容高，因此堆芯整體表現(xiàn)出較大的熱慣性，溫度變化滯后于功率變化。這些特性說明NUSTER-100的自穩(wěn)能力不強(qiáng)，這要求控制棒價(jià)值劃分需更加精細(xì)，步數(shù)更多，且對系統(tǒng)控制器要求更高。

圖24 +10 pcm反應(yīng)性引入工況下反應(yīng)堆的瞬態(tài)響應(yīng)

圖25為+10 pcm反應(yīng)性引入工況下堆芯各結(jié)構(gòu)平均溫度及系統(tǒng)電功率的瞬態(tài)響應(yīng)。由于各結(jié)構(gòu)間均為導(dǎo)熱關(guān)系，各層結(jié)構(gòu)互相之間的溫度響應(yīng)十分迅速，滯后性非常小。達(dá)到新穩(wěn)態(tài)后，電功率上升了3.65 kW左右，電功率的超調(diào)量約為1%。

圖25 +10 pcm反應(yīng)性引入工況下各結(jié)構(gòu)平均溫度和系統(tǒng)電功率的瞬態(tài)響應(yīng)

5 結(jié)論

本文提出了滿足無人水下潛航器安全可靠能源需求的海洋靜默式熱管反應(yīng)堆NUSTER-100概念設(shè)計(jì)方案，通過熱管反應(yīng)堆全系統(tǒng)數(shù)學(xué)物理模型的建立和自主知識產(chǎn)權(quán)瞬態(tài)分析程序HEART的開發(fā)，開展了熱管反應(yīng)堆穩(wěn)態(tài)和瞬態(tài)熱工水力特性研究，得到如下結(jié)論。

1)建立了熱管反應(yīng)堆全系統(tǒng)數(shù)學(xué)物理模型，開發(fā)了熱管反應(yīng)堆瞬態(tài)分析程序HEART，并通過了熱管實(shí)驗(yàn)、溫差發(fā)電實(shí)驗(yàn)等數(shù)據(jù)的驗(yàn)證與確認(rèn)，證明了模型和程序的正確性。

2)基于HEART開展的熱管反應(yīng)堆滿功率穩(wěn)態(tài)工況計(jì)算表明，熱管反應(yīng)堆的固態(tài)堆芯具有較好的溫度自展平能力，NUSTER-100熱電轉(zhuǎn)換效率達(dá)到10.88%。

3)冷態(tài)啟動瞬態(tài)工況研究發(fā)現(xiàn)，熱管反應(yīng)堆存在啟動速率限值，過快的升功率速度將導(dǎo)致反應(yīng)堆啟動失敗，因此提出了三段式啟動方案，可保證熱管反應(yīng)堆的平穩(wěn)啟動。

4)熱管反應(yīng)堆在反應(yīng)性引入瞬態(tài)工況時(shí)，堆芯整體表現(xiàn)出較大的熱慣性，溫度變化滯后于功率變化，核功率的超調(diào)量高于電功率的超調(diào)量，因此需要更加精細(xì)的控制棒價(jià)值劃分。