一種獲得偏離泡核沸騰型臨界熱流密度的方法
【專利摘要】一種獲得偏離泡核沸騰型臨界熱流密度的方法，步驟如下：1.設(shè)置一個初始熱流密度qm，為了獲得用于臨界熱流密度機(jī)理模型所需的局部熱工水力參數(shù)，基于均勻流模型建立管內(nèi)兩相熱工水力瞬態(tài)求解模型，基于有限差分法來求解兩相均勻流偏導(dǎo)，以及采用半隱式差分格式求解離散后的控制方程；2.當(dāng)步驟1計算所得局部熱工水力參數(shù)達(dá)到穩(wěn)態(tài)時，將步驟1獲得的局部熱工水力參數(shù)作為已知參數(shù)賦給偏離泡核沸騰型臨界熱流密度機(jī)理模型，求解臨界熱流密度qCHF；3.判斷qm與qCHF是否一致，若沒有收斂，則更新qm的值，重復(fù)步驟1，直至二者相同；本發(fā)明方法具有較好的拓展性，對數(shù)據(jù)依賴性低，能夠有效地預(yù)測臨界熱流密度，對反應(yīng)堆安全分析具有重要意義。
【專利說明】
一種獲得偏離泡核沸騰型臨界熱流密度的方法
技術(shù)領(lǐng)域
[0001] 本發(fā)明屬于反應(yīng)堆安全分析技術(shù)領(lǐng)域，具體涉及一種獲得偏離泡核沸騰型臨界熱 流密度的方法。 技術(shù)背景
[0002] 為保證堆芯流通通道的正常運(yùn)行和完整性，需要對管內(nèi)流動換熱特性進(jìn)行研究， 臨界熱流密度是反應(yīng)堆設(shè)計的主要熱工安全準(zhǔn)則之一。臨界熱流密度是核燃料元件表面發(fā) 生傳熱惡化時的熱通量，是冷卻劑流動沸騰機(jī)理發(fā)生轉(zhuǎn)變的結(jié)果。如果反應(yīng)堆燃料元件表 面發(fā)生沸騰臨界，將會導(dǎo)致燃料元件表面溫度過高從而造成加熱壁面燒毀，放射性物質(zhì)泄 露，造成嚴(yán)重的運(yùn)行事故，直接影響到反應(yīng)堆的安全運(yùn)行。因此準(zhǔn)確預(yù)測臨界熱流密度對反 應(yīng)堆的設(shè)計和安全分析、對反應(yīng)堆系統(tǒng)綜合性能的提高具有重要意義。
[0003] 目前，臨界熱流密度的預(yù)測方法較多采用適用于不同流體和不同結(jié)構(gòu)流道的各種 經(jīng)驗關(guān)系式以及臨界熱流密度查詢表。經(jīng)驗關(guān)系式是最早使用的預(yù)測和研究臨界熱流密度 的一種方法。至今，針對不同的熱工水力條件、流型、管道結(jié)構(gòu)尺寸等參數(shù)，國內(nèi)外研究者們 已經(jīng)提出數(shù)百種臨界熱流密度預(yù)測經(jīng)驗關(guān)系式。但大多數(shù)經(jīng)驗關(guān)系式是基于一個比較小的 數(shù)據(jù)庫開發(fā)得到，因此大多數(shù)經(jīng)驗關(guān)系式適用范圍有限，不能擴(kuò)展至適用范圍以外的工況。 而臨界熱流密度查詢表可在很寬的參數(shù)范圍直接查表并插值的方法得到相應(yīng)參數(shù)下的臨 界熱流密度，但其只適用于單一流體，且外推到查詢表提供的參數(shù)范圍之外(例如高流速和 高過冷度)時預(yù)測精度下降，較經(jīng)驗關(guān)系式需要花費(fèi)更多的計算時間。此外，數(shù)據(jù)表格的形 式使其不便于更新。
[0004] 基于此，有必要提供一種可避免使用大量的實(shí)驗數(shù)據(jù)擬合公式，降低對數(shù)據(jù)的依 賴性，并能夠適用于不同的流型和流體工質(zhì)的更有效率的預(yù)測臨界熱流密度的方法。
[0005] 發(fā)明說明
[0006] 為了解決上述問題，本發(fā)明提供一種獲得偏離泡核沸騰型臨界熱流密度的方法， 可避免使用大量的實(shí)驗數(shù)據(jù)擬合公式，降低對數(shù)據(jù)的依賴性，并能夠適用于不同的流型和 流體工質(zhì)。
[0007] 為了達(dá)到上述目的，本發(fā)明的采用如下技術(shù)方案：
[0008] -種獲得偏離泡核沸騰型臨界熱流密度的方法，包括如下步驟：
[0009] 步驟一:給定流道結(jié)構(gòu)尺寸以及邊界條件(進(jìn)口溫度、進(jìn)口流量和出口壓力），設(shè)置 一個初始熱流密度qm，采用均勻流模型模擬流道內(nèi)兩相瞬態(tài)流場，以此獲得流通通道中的 局部熱工水力參數(shù)，即壓力、流量及焓值，從而為偏離泡核沸騰型臨界熱流密度計算提供參 數(shù)；
[0010] 步驟二：當(dāng)步驟一中計算所得局部熱工水力參數(shù)達(dá)到穩(wěn)態(tài)時，將步驟一得到的出 口處局部熱工水力參數(shù)作為已知參數(shù)賦給偏離泡核沸騰型臨界熱流密度機(jī)理模型，采用該 機(jī)理模型計算得到臨界熱流密度qcHF;;
[0011]步驟三：比較步驟一中的初始熱流密度qm與步驟二計算所得臨界熱流密度qcHF，， 若qm大于qcHF，則減小qm，若qm小于qcHF，則增加 qm，采用變化后的熱流密度重復(fù)步驟一的過 程;若二者相同，則得到所求工況下的臨界熱流密度qcHF;
[0012] 步驟一所述的采用均勻流模型模擬流道內(nèi)兩相瞬態(tài)流場，采用有限差分法求解兩 相均勻流偏導(dǎo)，采用半隱式差分格式求解離散方程。
[0013] 步驟二所述的偏離泡核沸騰型臨界熱流密度機(jī)理模型為微液層蒸干模型，其原理 為:在壁面熱流密度逐漸緩慢增加過程中，流體達(dá)到兩相，壁面附近會出現(xiàn)一個拉長的汽 塊，它由多個小汽泡結(jié)合而成;長汽塊下裹覆著一層很薄的液體層(微液層），該長汽塊會沿 加熱壁面移動，此時長汽塊會阻礙主流液體對壁面的有效冷卻；當(dāng)長汽塊下的微液層蒸干 時，該點(diǎn)處的傳熱會劇烈惡化，從而導(dǎo)致壁溫急劇升高，進(jìn)而發(fā)生沸騰臨界;求解機(jī)理模型 的關(guān)鍵公式如下：
[0014]
[0015] 式中：
[0016] Pf--飽和液相密度/kg · m-3;
[0017] hfg--汽化潛熱/J .kg-1;
[0018] δ一一長汽塊下微液層厚度/m;
[0019] Ub--長汽塊移動速度/m · s-1;
[0020] Lb--長汽塊長度/m。
[0021] 從式（1)中可知，微液層厚度(δ)、汽塊移動速度(Ub)和汽塊長度(Lb)是求解的關(guān)鍵 參數(shù)；
[0022 ]其中汽塊移動速度Ub通過求解汽塊的軸向受力平衡來確定，軸向方向上汽塊的受 力包括浮力Fb和拖曳力Fd，公式如下：
[0023] Fb-Fd = O 公式（2)
[0024]微液層厚度δ通過求解汽塊的徑向受力平衡來確定，徑向方向上的汽塊受力平衡 包括側(cè)面提升力Fr、蒸發(fā)力Fe3以及壁面潤滑力Fw:
[0025] -Fr+Fe+Fw=0 公式（3)
[0026]汽塊長度Lb采用Helmholtz臨界波長，其公式如下所示：
[0027]
[0028] 式中：
[0029] P1--液相密度/kg · m-3;
[0030] pg一一汽相密度/kg· m-3;
[0031 ] σ--表面張力/N · m工。
[0032] 本發(fā)明和現(xiàn)有技術(shù)相比，具有如下優(yōu)點(diǎn)：
[0033] 1.采用均勻流模型模擬了兩相瞬態(tài)流場，可較精確地計算模擬出臨界沸騰發(fā)生點(diǎn) 處的熱工水力參數(shù)。
[0034] 2.從微觀上對汽泡受力平衡和運(yùn)動進(jìn)行研究，得到了較為精細(xì)的機(jī)理模型，同時 也便于后續(xù)分析運(yùn)動條件對模型的影響因素。
【附圖說明】
[0035] 圖1為本發(fā)明方法的流程圖。
[0036] 圖2為一維空間網(wǎng)格。
[0037]圖3為汽塊受力的示意圖。
[0038]圖4為機(jī)理模型計算流程圖。
【具體實(shí)施方式】
[0039]下面結(jié)合附圖和【具體實(shí)施方式】對本發(fā)明結(jié)構(gòu)進(jìn)行詳細(xì)說明。
[0040]如圖1所示，本發(fā)明一種獲得偏離泡核沸騰型臨界熱流密度的方法，具體方法如 下：
[0041 ]步驟一:給定流道結(jié)構(gòu)尺寸以及邊界條件(進(jìn)口溫度、進(jìn)口流量和出口壓力），設(shè)置 一個初始熱流密度qm，采用均勻流模型模擬流道內(nèi)兩相瞬態(tài)流場，以此獲得流通通道中的 局部熱工水力參數(shù)，即壓力、流量及焓值，從而為偏離泡核沸騰型臨界熱流密度計算提供參 數(shù)；
[0042] 步驟二：當(dāng)步驟一中計算所得局部熱工水力參數(shù)達(dá)到穩(wěn)態(tài)時，將步驟一得到的出 口處局部熱工水力參數(shù)作為已知參數(shù)賦給偏離泡核沸騰型臨界熱流密度機(jī)理模型模塊，計 算得到臨界熱流密度qcHF;
[0043]步驟三：比較步驟一中的初始熱流密度qm與步驟二計算所得臨界熱流密度qcHF，若 qm大于qcHF，則減小qm，若qm小于qcHF，則增加 qm，采用變化后的熱流密度重復(fù)步驟一的過程； 若二者相同，則得到所求工況下的臨界熱流密度qcHF。
[0044] 下面進(jìn)行詳細(xì)說明：
[0045] 步驟一中：兩相瞬態(tài)流場采用均勻流模型模擬求解，下面給出具體的模型：
[0046] 由于所研究的臨界熱流密度特性只需要提供外部的熱工水力參數(shù)，對流場模擬的 要求不高，因此采用改進(jìn)的均勻流模型模擬流場。均勻流模型的基本假設(shè)為：
[0047]①兩相具有相等的線速度，即Ug = Ul = UH;
[0048]②兩相之間處于熱力平衡；
[0049]③使用合理確定的單相摩擦因數(shù)。
[0050]運(yùn)用上述假定于兩相流動基本方程，并與單相流體守恒方程類比，定義合適的混 合物物性和摩擦因數(shù)計算式，可以得到以下守恒方程。
[0051 ] 1)質(zhì)量守恒方程：
[0052]
公式（5)
[0053] 式中：
[0054] p一一兩相混合物的密度/kg · Hf3;
[0055] u一一兩相混合物的流動速度/m · ?Γ1;
[0056] t--時間/s;
[0057] z--長度/m。
[0058] 考慮兩相流體的可壓縮性，為了求解方程組，增加如下物性方程：
[0059] p = p(p,h)公式（6)
[0060]由式(6)求導(dǎo)可得：
[0061]
[0062]
[0063]
[0064]
[0065]
[0066]
[0067]
[0068]
[0069]
[0070]
[0071]
[0072] 式中：
[0073] Prh 加熱周長/m;
[0074] A--流通面積/m2;
[0075] qv--體積釋熱率/W · m 3;
[0076] ①若無內(nèi)熱源，即流體內(nèi)部沒有釋熱。上式右側(cè)的qv為零。
[0077] ②忽略動能變化引起的能量變化，則上式左偵
[0078]因此，能量方程可簡化為
[0079]
[0080]
[0081]
[0082]
[0083]
[0084]步驟一中：上述守恒方程，采用有限差分法求解偏導(dǎo)，采用半隱式差分格式求解離 散方程。
[0085]在空間上，控制體劃分采用交錯網(wǎng)格。交錯網(wǎng)格中，熱力學(xué)變量(壓力p，焓值h，密 度P，空泡份額α等)存儲在控制體中心（圖2中虛線位置），而流體的速度存放在控制體邊界 處(圖2中實(shí)線位置）。
[0086]在時間上采用半隱式差分格式進(jìn)行離散。為了加快計算速度，除了連續(xù)方程和能 量方程中對流項的速度，以及動量方程中的壓力梯度項采用隱式格式以外，其他項均采用 顯示格式。
[0087] 各守恒方程的離散形式如下：
[0088] 席醫(yī)Φ栢卞'fS.
[0089]
[0090]
[0091]
[0092]
[0093]
[0094]
[0095] 時間步長的選取可以違反(u+a) Δ t/ Δ Z<1的波速限制，但是要滿足u Δ t/ Δ Z<1 的流速限制。
[0096] 能量守恒方程：
[0097]
[0098]
[0099]
[0100] 其中Prh 加熱周長/m。
[0101] 交錯網(wǎng)格以及有限差分法的優(yōu)點(diǎn)是壓力梯度項與速度項自然的耦合在一起，并且 如果滿足Courant法則u Δ t/ Δ Z<1則形成一個對角優(yōu)勢矩陣，可保證求解算法的穩(wěn)定性。 所有網(wǎng)格的方程最后可以組成一個壓力的三對角矩陣，它可以利用直接消去法求解，也可 以利用Gauss-Seidel迭代法求解?？傮w計算思路是:①用動量方程求解出速度;②利用上一 時層求解壓力方程;③用解出的壓力回代求解速度和焓值。
[0102] 步驟二中：微液層蒸干機(jī)理的具體模型如下：
[0103] 微液層蒸干模型假設(shè)加熱壁面附近產(chǎn)生的小汽泡結(jié)合形成大汽塊，在汽塊下存在 非常薄的液相層，稱為微液層。汽塊移動過程中，當(dāng)汽塊下的液體全部蒸發(fā)燒干時，該點(diǎn)處 的加熱壁面被單相蒸汽覆蓋從而導(dǎo)致傳熱惡化，進(jìn)而導(dǎo)致沸騰臨界發(fā)生。求解機(jī)理模型的 關(guān)鍵公式如下：
[0104] qcHF = 5pfhfgUb/Lb 公式（19)
[0105] 式中：
[0106] pf--飽和液相密度/kg · m-3;
[0107] hfg一一汽化潛熱/J .kg-1;
[0108] δ一一汽塊下微液層厚度/m;
[0109] Ub--汽塊移動速度/m · s-1;
[0110] Lb--汽塊長度/m。
[0111] 從（14)式中可知，微液層厚度(δ)、汽塊移動速度(Ub)和汽塊長度(Lb)是求解的關(guān) 鍵參數(shù)。
[0112]如圖3所示，上述三個關(guān)鍵參數(shù)可以通過分析汽塊受力獲得，其中汽塊移動速度Ub 通過求解汽塊的軸向受力平衡來確定。軸向方向上汽塊的受力包括浮力Fb和拖曳力Fd，公式 如下：
[0113]
[0114]
[0115]
[0116]
[0117] 式中：
[0118] Db--汽塊當(dāng)量直徑/m;
[0119] P1--液相的密度/kg · m-3;
[0120] pg一一汽相的密度/kg· m-3;
[0121] g一一重力加速度/miT2;
[0122] t--時間/s;
[0123] Cd--拖拽系數(shù)；
[0124] Ubi一一汽塊中心線所處的徑向位置處的主流速度/ms'
[0125] 聯(lián)立式(20)-(22)，汽塊速度Ub可表達(dá)如下：
[0126]
[0127]
[0128]
[0129]
[0130]
[0131] 式中：
[0132] σ--表面張力/N · m S
[0133] D--管道直徑/m;
[0134] Tw--壁面剪切力/kg · m-1 · s-2。
[0135 ] 采用Karman速度分布方程確定當(dāng)?shù)亓魉伲?br>[0136]
[0138] 近壁面區(qū)的汽塊一般處在緩沖區(qū)的范圍，即對應(yīng)于Karman速度分布方程的第二個 關(guān)系式。所以為簡化計算，可選用第二個關(guān)系式計算當(dāng)?shù)亓魉俜植?，該關(guān)系式可以改寫為：
[0137]
[0139]
[0140]
[0141]
[0142]
[0143]
[0144] 微液層厚度δ通過求解汽塊的徑向受力平衡來確定。徑向方向上的汽塊受力包括 平衡側(cè)面提升力Fr、蒸發(fā)力Fe3以及壁面潤滑力F w:
[0145] -Fr+Fe+Fw=0 公式（30)
[0146] 蒸發(fā)力是由微液層中液相蒸發(fā)產(chǎn)生的小汽泡進(jìn)入汽塊，對汽塊沖擊而產(chǎn)生的，這 個力將推動汽塊向遠(yuǎn)離加熱壁面。側(cè)面提升力是由汽塊和主流區(qū)存在速度差異，使汽塊在 運(yùn)動過程中發(fā)生旋轉(zhuǎn)而產(chǎn)生，提升力的大小還與汽塊在徑向方向上的速度梯度有關(guān)，它會 推動汽塊靠近加熱壁面。
[0147] 沸騰臨界將要發(fā)生時，由微液層蒸發(fā)而產(chǎn)生的蒸發(fā)力為：
[0148]
公Λ (31)[0149] 垂直管道內(nèi)兩相湍流流動中的側(cè)面提升力的表達(dá)式：
[0150；
[0151；
[0152；
[0153；
[0154] 其中式(33)中，α為臨界沸騰發(fā)生點(diǎn)處的空泡份額，式(29)中的速度由式(22)計算 獲得。
[0155] 因此，側(cè)面提升力可表達(dá)如下：
[0156]
[0157] 此外還需要考慮壁面潤滑力Fw，壁面潤滑力將非常靠近壁面處的汽泡推向中心 區(qū)，而對離壁面距離較遠(yuǎn)的汽泡作用力則非常小，其表達(dá)式如下：
[0158]
[0159]
[0160]
[0161] 上式中，&和(：2為經(jīng)驗常數(shù)，yw表示離加熱壁面的垂直距離，式(36)在y w〈5Db范圍內(nèi) 有效。
[0162] 聯(lián)立式（30)，（31)，（32)和(35)，微液層厚度δ可表達(dá)如下：
[0163]
[0164] 求出δ以后，可以通過式(28)計算得到Ubl，然后分別用式(23)和(29)計算得到U b和 U，再將新的Ub和Lb代入式(38)計算得到新的δ，通過迭代計算直到以上幾個參數(shù)值達(dá)到收 斂。最后，將迭代計算得到的S，U b和Lb代入式（19)來計算臨界熱流密度值。機(jī)理模型的流程 圖如圖4所示。
【主權(quán)項】
1. 一種獲得偏離泡核沸騰型臨界熱流密度的方法，其特征在于:包括如下步驟： 步驟一:給定流道結(jié)構(gòu)尺寸以及邊界條件即進(jìn)口溫度、進(jìn)口流量和出口壓力，設(shè)置一個 初始熱流密度qm，采用均勻流模型模擬流道內(nèi)兩相瞬態(tài)流場，以此獲得流通通道中的局部 熱工水力參數(shù)，即壓力、流量及焓值，從而為偏離泡核沸騰型臨界熱流密度計算提供參數(shù)； 步驟二：當(dāng)步驟一中計算所得局部熱工水力參數(shù)達(dá)到穩(wěn)態(tài)時，將步驟一得到的出口處 局部熱工水力參數(shù)作為已知參數(shù)賦給偏離泡核沸騰型臨界熱流密度機(jī)理模型，采用該機(jī)理 模型計算得到臨界熱流密度qcHF; 步驟三:比較步驟一中的初始熱流密度qm與步驟二計算所得臨界熱流密度qcHF，，若qm大 于qCHF，則減小qm，若qm小于qCHF，則增加 qm，采用變化后的熱流密度重復(fù)步驟一的過程;若二 者相同，則得到所求工況下的臨界熱流密度qcHF。2. 根據(jù)權(quán)利要求1所述的一種獲得偏離泡核沸騰型臨界熱流密度的方法，其特征在于： 步驟一所述的采用均勻流模型模擬流道內(nèi)兩相瞬態(tài)流場，采用有限差分法求解兩相均勻流 偏導(dǎo)，采用半隱式差分格式求解離散方程。3. 根據(jù)權(quán)利要求1所述的一種獲得偏離泡核沸騰型臨界熱流密度的方法，其特征在于： 步驟二所述的偏離泡核沸騰型臨界熱流密度機(jī)理模型為微液層蒸干模型，其原理為:在壁 面熱流密度逐漸緩慢增加過程中，流體達(dá)到兩相，壁面附近會出現(xiàn)一個拉長的汽塊，它由多 個小汽泡結(jié)合而成；長汽塊下裹覆著一層很薄的液體層即微液層，該長汽塊會沿加熱壁面 移動，此時長汽塊會阻礙主流液體對壁面的有效冷卻；當(dāng)長汽塊下的微液層蒸干時，該點(diǎn)處 的傳熱會劇烈惡化，從而導(dǎo)致壁溫急劇升高，進(jìn)而發(fā)生沸騰臨界;求解機(jī)理模型的關(guān)鍵公式 如下： qcHF = 8pfhfgUb/Lb 公式（1) 式中： Pf--飽和液相密度/kg ? nf3; hfg一一汽化潛熱/J ? kg' 8一一長汽塊下微液層厚度/m; Ub--長汽塊移動速度/m ? s S Lb--長汽塊長度/m。 從式(1)中可知，微液層厚度S、汽塊移動速度Ub和汽塊長度Lb是求解的關(guān)鍵參數(shù)； 其中汽塊移動速度Ub通過求解汽塊的軸向受力平衡來確定，軸向方向上汽塊的受力包 括浮力Fb和拖曳力Fd，公式如下： Fb-Fd = 0 公式(2) 微液層厚度S通過求解汽塊的徑向受力平衡來確定，徑向方向上的汽塊受力平衡包括 側(cè)面提升力Fr、蒸發(fā)力以及壁面潤滑力Fw: -Fr+Fe+Fw=0 公式(3) 汽塊長度Lb采用Helmholtz臨界波長，其公式如下所示：Pi一一液相密度/kg ? m-3; Pg一一汽相密度/kg ? nf3; 8--表面張力/N ? nf1。
【文檔編號】G06F17/50GK106055850SQ201610566357
【公開日】2016年10月26日
【申請日】2016年7月18日
【發(fā)明人】蘇光輝, 黃思洋, 張大林, 田文喜, 秋穗正
【申請人】西安交通大學(xué)