侯 東,許志紅,林 萌,劉鵬飛,楊燕華
(上海交通大學核科學與工程學院,上海 200240)
超臨界水堆(SCWR)是第四代核反應堆中唯一的壓水堆,具有熱效率高、系統(tǒng)簡單、安裝尺寸小和經濟性好等眾多優(yōu)點。但是,由于堆芯進出口之間超臨界水具有較大的焓差致使其密度發(fā)生劇烈變化,并容易引起系統(tǒng)振蕩而導致發(fā)生不穩(wěn)定現象。因此,加熱通道內超臨界水的流動穩(wěn)定性問題引起了國內外學者的廣泛關注[1-8]。
目前國際上主要通過將一維單相流動偏微分控制方程直接進行線性化處理,然后應用頻域法在給定的設計參數和工況下對SCWR堆芯進行穩(wěn)定性分析[2],通過對系統(tǒng)衰減率[2]的計算來考察系統(tǒng)是否穩(wěn)定。這種方法不具有普適的通用性,只能對特定的系統(tǒng)參數給出穩(wěn)定性判定。而本文則主要是從超臨界壓力下水的熱工水力學特性本身出發(fā),對垂直加熱通道內超臨界水流動穩(wěn)定性機理進行研究,從而得出廣泛適用的穩(wěn)定性判斷準則,為后續(xù)相關的SCWR設計提供指導。本文將無量綱分析法和頻域分析法相結合,定義了影響系統(tǒng)穩(wěn)定性的關鍵參數,并給出了各種不同參數所對應的穩(wěn)定性邊界。同時,所提出的無量綱頻域分析方法也可為廣大科研工作者提供參考。
國際上現在主要存在兩大類超臨界流動分析方法:一類是將超臨界流體流動看作是類似汽液兩相流,具體的兩相分界點在不同的文獻中有不同的定義[2,4];另一類是將超臨界流動看作是單相流動,但在流動過程中物性變化比較劇烈。目前國際上沒有統(tǒng)一的將超臨界流體劃分成兩相的標準,同時兩相流模型中通常需要對流體物性進行簡化處理[2],使得分析結論具有很大的不確定度。雖然采用單相流模型不能很好的處理系統(tǒng)壓力跨流體臨界壓力時的相變問題,但是本文的重點是研究大于臨界壓力的系統(tǒng)穩(wěn)定性,因此可以不用考慮系統(tǒng)壓力降低到臨界壓力下時的相變情況。同時,相對于兩相流模型,單相流模型具有形式簡單的優(yōu)點,不需要對流體物性進行特殊簡化,能更準確的描述超臨界流體的流動。因此本文采用單相流動模型來分析超臨界流體的流動。
現有的公開文獻中大都只對無量綱分析方法進行基本闡述,而未給出具體的實現方法。本文將給出無量綱分析與頻域法分析相結合的分析方法的具體實施步驟,以為廣大科研工作者提供借鑒。
進行無量綱分析的第一步是將原始偏微分控制方程進行無量綱化。管內一維單相流動控制方程(質量守恒方程、動量守恒方程、能量守恒方程)分別為:
式中:ql是平均線功率;δz0和 δzL分別在z=0和z=L的時候為1,z為其他值時為0;Dh為管道水力直徑;f為管道壁面摩擦因數;Kin和Kout分別為管道進出口阻力因數;A為管道流通面積;t表示時間;z表示空間位置;其他參數與常規(guī)熱工水力參數意義相同。
定義如下無量綱數[4]:
式中 :ρPC、hPC、Cp,PC和 βPC分別表示給定系統(tǒng)壓力下的擬臨界溫度點的流體密度、焓、定壓比熱容和熱膨脹系數;μin表示入口流體流速;Fr為Frout數;其他參數與常規(guī)熱工水力參數意義相同。
將各無量綱數代入方程(1)~(3)中,可得無量綱一維流動控制方程。無量綱質量守恒方程:
無量綱動量守恒方程:
無量綱能量守恒方程:
其中無量綱參數NSPC和NTPC是兩個重要的確定系統(tǒng)穩(wěn)定性的參數。從它們的定義可知NSPC主要表征系統(tǒng)入口流體的無量綱過冷度,即系統(tǒng)的入口流體狀態(tài);而NTPC主要表征了流體經過系統(tǒng)加熱后的體積膨脹率,因此某些文獻中[2]將其稱為流體膨脹數。系統(tǒng)穩(wěn)定性邊界是由這兩個參數組成的點集構成[5]。
經上述無量綱化后,無量綱密度可近似成無量綱焓的單值函數,即無量綱狀態(tài)方程可表示為[4]:
以壓力分別為 24.5 MPa、25.0 MPa和30.0 MPa時的無量綱密度和無量綱焓為例,如圖1所示,三條曲線完全重合。從國外學者的研究文獻中可以發(fā)現,包括其他非水工質(如CO2等)在內的各種流體的無量綱密度和無量綱焓的關系曲線也是基本重合的[5]。因此,上述無量綱化處理具有廣泛的適用性。
將方程(4)代入到方程(6)中可得:
圖1 無量綱密度隨無量綱焓的變化Fig.1 Dimensionless density vs.dimensionless enthalpy
方程(5)、(6)或(7)、(9)和狀態(tài)方程(8)一起構成無量綱控制方程。
為消除空間微分算符,需要對無量綱控制方程進行空間離散,本文在軸向上進行均勻劃分網格,采用有限控制體積法進行離散。其中質量與能量方程的控制體是重合的,而動量方程是按照交錯網格的方式進行離散化的,其控制體相對質量和能量方程的控制體偏移半個控制體積??臻g離散方程如下:
將狀態(tài)方程代入到質量方程中,可得
上述擾動方程組中各擾動系數由穩(wěn)態(tài)計算值所確定。
將線性擾動方程(13)~(15)表示成如下矩陣形式為:
式中:
A和B分別為線性擾動方程的時間導數項系數矩陣和空間離散項系數矩陣。根據頻域法的理論[2,6-7],對式(16)時間導數項進行 Lap lace變換可得:
系統(tǒng)的特征方程為:
該方程實部最大的根為系統(tǒng)的主特征根λRe+iλIm,根據其實部 λRe便可確定系統(tǒng)的穩(wěn)定狀態(tài)。根據頻域法理論可知,當λRe>0時系統(tǒng)是非穩(wěn)定的,當λRe<0時系統(tǒng)是穩(wěn)定的,而當λRe=0時系統(tǒng)處于臨界穩(wěn)定狀態(tài)。
從系統(tǒng)特征方程(18)可以看出,在計算矩陣的行列式時矩陣中含有非數值元素s,因此不能用常規(guī)的數值求解方法來求解式(18)。本文是利用M ATLAB軟件提供的符號運算功能[10]來進行特征方程的建立和求解的。
從上述系統(tǒng)特征方程的推導可知,確定系統(tǒng)的穩(wěn)定邊界前必須先確定系統(tǒng)的穩(wěn)態(tài)參數,然后才能求出系統(tǒng)主特征根并確定系統(tǒng)的穩(wěn)定性[6-8]。
在穩(wěn)態(tài)情況下,所有無量綱方程的時間導數項均為0,根據空間離散方程(10)~(12)有:
式中:G*為NSPC所確定,各節(jié)塊h*也能根據能量方程(21)確定,進而ρ*也可以由狀態(tài)方程(8)確定。根據動量方程(20)可求得進出口壓差為Fr的函數,可表示成如下形式:
式中:a和b是根據方程(20)計算得到的常數。
由于本文考查的是流動穩(wěn)定性中較為重要的密度波穩(wěn)定性,所計算的問題中進出口壓差是給定的常數,即 ΔP*為常數,因此可根據方程(22)可以確定Fr值,進而可以求出系統(tǒng)特征方程的系數矩陣。
當穩(wěn)態(tài)參數分布確定后,擾動方程的系數矩陣便可以確定,從而系統(tǒng)的特征方程也可確定。在給定NSPC時,通過改變NTPC的值使主特征根的最大實部為 0,由此得到的點(NSPC,NTPC)為一個穩(wěn)定邊界點。所考察范圍內的所有邊界點一起組成一條穩(wěn)定邊界線,由此可以得到系統(tǒng)的穩(wěn)定邊界。
本文所計算問題的系統(tǒng)壓力為25MPa,因為目前國內外所設計的SCWR系統(tǒng)壓力均為25 MPa。其他相關基本參數如表1所示。
其中,重力加速度g為負數是因為本文主要考慮向下流動加熱管道內的流動穩(wěn)定性,當流體向上流動時g應為正數。
本文通過對不同的系統(tǒng)進口阻力因數Kin、出口阻力因數Kout、摩擦因數 Λ、進出口壓差ΔP的穩(wěn)定性邊界計算來觀察各參數對超臨界流體流動的穩(wěn)定性影響。上海交通大學設計的混合能譜SCWR堆芯內流體是在下降加熱通道中被加熱到擬臨界溫度以上的[9],因此本文主要考慮垂直向下流動的情況,但也同時與流動方向向上時的系統(tǒng)穩(wěn)定性邊界進行了對比。
表1 穩(wěn)定性分析系統(tǒng)參數Table1 Stability analysis system parameters
計算了入口阻力因數分別為0和20時的系統(tǒng)穩(wěn)定性邊界,不同入口阻力系數下的穩(wěn)定性邊界如圖2所示。高的入口阻力因數有利于系統(tǒng)的穩(wěn)定。
圖2 不同入口阻力因數下的穩(wěn)定性邊界Fig.2 The stability boundary at different inlet orifice coefficient
計算了出口阻力因數分別為1和20時的系統(tǒng)穩(wěn)定性邊界,不同出口邊界阻力因數下的穩(wěn)定性邊界如圖3所示。高的出口阻力因數不利于系統(tǒng)的穩(wěn)定。
圖3 不同出口阻力因數下的穩(wěn)定性邊界Fig.3 The stability boundary at different outlet orifice coefficient
計算了無量綱摩擦因數分別為5和22時的系統(tǒng)穩(wěn)定性邊界,不同摩擦因數下的穩(wěn)定性邊界如圖4所示。高的摩擦因數不利于系統(tǒng)的穩(wěn)定。
圖4 不同摩擦因數下的穩(wěn)定性邊界Fig.4 The stability boundary at different friction coefficient
計算了系統(tǒng)壓差分別為 0.15 MPa和0.25 MPa時的系統(tǒng)穩(wěn)定性邊界,不同進出口壓差下的穩(wěn)定性邊界如圖5所示。系統(tǒng)進出口壓差對系統(tǒng)的穩(wěn)定性影響較小。
圖5 不同系統(tǒng)壓差下的穩(wěn)定性邊界Fig.5 The stability boundary at different system pressure drop
計算了垂直向上和垂直向下流動的穩(wěn)定性邊界,如圖6所示。流體向上流動比向下流動更有利于系統(tǒng)的穩(wěn)定。
圖6 流動穩(wěn)定性邊界Fig.6 Flow stability boundary
以上穩(wěn)定性邊界計算分析可以對系統(tǒng)的穩(wěn)定性進行初步判斷,為SCWR的堆芯和系統(tǒng)設計提供指導。
(1)給出了將無量綱分析與頻域法分析相結合來確定系統(tǒng)穩(wěn)定性邊界的詳細描述,包括方法論和具體的實施步驟。
(2)通過對加熱通道內超臨界流體流動方程進行無量綱化、線性化,以及Lap lace變換處理,定義了超臨界流動的兩個重要無量綱參數NSPC和NTPC,并以此畫出了系統(tǒng)穩(wěn)定邊界圖。
(3)對系統(tǒng)入口阻力因數、出口阻力因數、摩擦因數、進出口壓降,以及流動方向等進行了穩(wěn)定敏感性分析,為SCWR的堆芯和系統(tǒng)設計提供了指導。
本文下一步的工作將進行SCWR堆芯穩(wěn)定性分析,并考慮功率反饋下的堆芯響應,以期得到更系統(tǒng)化的穩(wěn)定性分析結論。
[1] 薛愛軍,程旭.超臨界水熱力系統(tǒng)的穩(wěn)定性的簡化模型分析[J].核動力工程,2009,30(5):35-38.
[2] Zhao J Y.Stability Analysis of Supercritical Water cooled Reactors[D].Departmen t of Nuclear Science and Engineering,M IT,Camb ridge,Massachusetts,2005.
[3] Lahey R T,Moody F J.The thermal-hydraulics of a boiling w ater reactor[J].American Nuclea Society,1993.
[4] Amb rosiniW.On the Analogies in the Dynamic Behavior of H eated Channels with Boiling and Supercritical Fluids[J].Nu clear Engineering and Design,2007,237:1 164-1 174.
[5] Amb rosiniW,Sharabi M B.Assessment of stability maps fo r heated channelds with supercritical fluids versus the p redictions of a system code[J].Nuclear Engineering and Design,2007,39:627-635.
[6] Chatoorgoon V,Voodi A,Fraaser D.The Stability Boundary fo r Supercritical Flow in Natu ral Convection Loops Part 1:H 20 Studies[J].Nuc lear Engineering and Design,2005,235:2 570-2 580.
[7] A rgonne National Laboratory.Initial Imp lementation of Multi-Channel Therm al-H ydrau lics Capability in Frequen cy Domain SCWRStability Analy sis Code SCWRSA[Z].2005.
[8] Cheng X,Yang Y H.A Point-Hydraulic Model for Flow Stability Analy sis[J].Nuc lear Engineering and Design.2007,238:188-199.
[9] Liu X J,Cheng X.Thermal-hydrau lic and neu tronphysical characteristics of a new SCWR fuel assem bly[J].Annalsof Nu clear Energy,2009,36:28-36.
[10]The MathW orks Inc.M ATLAB Symbolic Math Toolbox 5.4[M].USA:The MathW orks Inc.,2000.