侯曉凡，孫中寧

（哈爾濱工程大學(xué)核安全與仿真技術(shù)國(guó)防重點(diǎn)學(xué)科實(shí)驗(yàn)室，黑龍江哈爾濱150001）

開(kāi)式自然循環(huán)閃蒸不穩(wěn)定的線性均相流模型

侯曉凡，孫中寧

（哈爾濱工程大學(xué)核安全與仿真技術(shù)國(guó)防重點(diǎn)學(xué)科實(shí)驗(yàn)室，黑龍江哈爾濱150001）

針對(duì)開(kāi)式自然循環(huán)閃蒸不穩(wěn)定邊界的判定問(wèn)題，本文運(yùn)用狀態(tài)空間方法對(duì)其進(jìn)行研究，提出了描述該過(guò)程的線性均相流模型，導(dǎo)出了判斷系統(tǒng)穩(wěn)定性的狀態(tài)空間表達(dá)式，并求解得到了閃蒸不穩(wěn)定的邊界。模型計(jì)算結(jié)果和實(shí)驗(yàn)結(jié)果符合較好，僅在高加熱功率工況下存在一定偏差，而形成偏差的原因是加熱段內(nèi)沸騰噴發(fā)過(guò)程對(duì)閃蒸過(guò)程產(chǎn)生強(qiáng)烈的干擾，使其無(wú)法形成穩(wěn)定的閃蒸流動(dòng)。同時(shí)實(shí)驗(yàn)中通過(guò)注氣誘導(dǎo)法排除了對(duì)閃蒸過(guò)程的影響，得到新的閃蒸不穩(wěn)定邊界，與模型計(jì)算得到的下邊界符合很好，證明了本模型結(jié)果具有較高的準(zhǔn)確性。

開(kāi)式自然循環(huán)；閃蒸不穩(wěn)定；線性均相流模型；狀態(tài)空間方程；注氣誘導(dǎo)法；穩(wěn)定閃蒸

流動(dòng)不穩(wěn)定是制約自然循環(huán)系統(tǒng)安全穩(wěn)定運(yùn)行的重要因素［1］，特別是開(kāi)式自然循環(huán)系統(tǒng)，由于其系統(tǒng)壓力較低，易于發(fā)生沸騰、流動(dòng)閃蒸等劇烈的相變，同時(shí)缺乏能夠容納系統(tǒng)壓力波動(dòng)，抑制空泡份額變化及流量振蕩的封閉氣空間，因此開(kāi)式自然循環(huán)系統(tǒng)更容易發(fā)生流動(dòng)不穩(wěn)定現(xiàn)象［2］。

在開(kāi)式自然循環(huán)流動(dòng)不穩(wěn)定中，閃蒸不穩(wěn)定是重要的不穩(wěn)定類型之一。閃蒸不穩(wěn)定發(fā)生時(shí)，上升段內(nèi)出現(xiàn)間歇性的閃蒸流動(dòng)或者閃蒸起始點(diǎn)沿上升通道周期性地上下移動(dòng)，從而導(dǎo)致自然循環(huán)流動(dòng)出現(xiàn)振蕩。目前，針對(duì)開(kāi)式自然循環(huán)閃蒸不穩(wěn)定的研究多為實(shí)驗(yàn)研究［3-6］，理論分析研究也集中于采用數(shù)值解析方法進(jìn)行求解［7-8］，而采用近似分析方法對(duì)閃蒸不穩(wěn)定進(jìn)行研究的文獻(xiàn)很少［9］。本文提出了描述開(kāi)式自然循環(huán)系統(tǒng)閃蒸流動(dòng)不穩(wěn)定的線性均相流模型，并基于系統(tǒng)控制原理中的狀態(tài)空間分析方法對(duì)閃蒸不穩(wěn)定邊界進(jìn)行分析。

1　數(shù)學(xué)模型建立

1.1開(kāi)式自然循環(huán)閃蒸流動(dòng)物理過(guò)程及基本假設(shè)

圖1給出了開(kāi)式自然循環(huán)閃蒸流動(dòng)時(shí)，加熱段和上升段內(nèi)的溫度壓力的變化情況。加熱段出口流體處于過(guò)冷狀態(tài)，隨著流體向上流動(dòng)，其承受的重位壓頭和相應(yīng)的飽和溫度均不斷降低，導(dǎo)致流體運(yùn)動(dòng)至上升段某一位置時(shí)達(dá)到當(dāng)?shù)仫柡蜏囟榷焖倨纬善簝上嗔鲃?dòng)，該位置稱為閃蒸起始點(diǎn)。

系統(tǒng)出現(xiàn)穩(wěn)定閃蒸流動(dòng)時(shí)，閃蒸起始點(diǎn)固定在上升段的某一位置上，而系統(tǒng)呈現(xiàn)閃蒸不穩(wěn)定特性時(shí)，閃蒸起始點(diǎn)間歇性地消失或沿上升段周期性上下移動(dòng)，當(dāng)閃蒸起始點(diǎn)消失（或下移）時(shí)，上升段內(nèi)空泡份額降低，自然循環(huán)驅(qū)動(dòng)降低，使得自然循環(huán)流量降低，上升段內(nèi)流體溫度增加，進(jìn)而導(dǎo)致閃蒸起始點(diǎn)出現(xiàn)（或上移）。隨后上升段內(nèi)空泡份額升高，循環(huán)流量上升，上升段內(nèi)流體溫度降低，又導(dǎo)致閃蒸起始點(diǎn)再次消失（或下移）。如此反復(fù)形成了自持性的閃蒸流動(dòng)不穩(wěn)定。

圖1　上升段內(nèi)閃蒸流動(dòng)示意圖Fig.1 Sketch of flashing flow in upward channel

為了建立自然循環(huán)系統(tǒng)閃蒸不穩(wěn)定的線性均相流模型，本文針對(duì)上述物理過(guò)程提出如下假設(shè):

1）流體為一維流動(dòng)，流型的影響可忽略不計(jì)；

2）上升段內(nèi)流動(dòng)過(guò)程（包括閃蒸過(guò)程）為熱力學(xué)平衡過(guò)程，同時(shí)忽略加熱段內(nèi)過(guò)冷沸騰的影響；

3）考慮到系統(tǒng)內(nèi)溫度和壓力的變化范圍較小，故忽略溫度和壓力變化對(duì)于液相密度ρl和氣相密度ρg的影響。

1.2基本控制方程及邊界條件

針對(duì)上升段內(nèi)的閃蒸流動(dòng)過(guò)程，列出均相流的質(zhì)量守恒方程和動(dòng)量守恒方程:

式中:t為時(shí)間，ρ為密度，α為空泡份額，G為質(zhì)量流速，z為距加熱段入口距離，P為壓力，D為當(dāng)量直徑，?f為分液相折算系數(shù)，ξ為摩擦阻力系數(shù)，下標(biāo)g表示氣相，下標(biāo)l表示液相。同時(shí)列出能量守恒方程時(shí)，忽略了外界力做功及動(dòng)能變化引起的能量變化:

式中:q為線功率密度，A為流道橫截面積，hl和hg分別為液相焓值和氣相焓值，x為干度。

另外，根據(jù)漂移流模型得到上升段兩相區(qū)內(nèi)流體的質(zhì)量含氣率x和空泡份額α之間的關(guān)系，列出輔助方程:

式中:C0為分布參數(shù)；υgj為氣相漂移速度，其具體數(shù)值的選取參照Z(yǔ)uber&Findlay公式［10］。

最后列出了開(kāi)式自然循環(huán)上升段內(nèi)閃蒸流動(dòng)的邊界方程，即

加熱段入口位置:

上升段入口位置:

上升段出口位置:

式中:下標(biāo)a為上升段出口位置，b為上升段入口位置，c表示加熱段區(qū)域，e為外界環(huán)境壓力，in和ex分別為入口和出口位置，f表示飽和狀態(tài)的參數(shù)值，xc，in為相對(duì)于上升段出口處飽和狀態(tài)的加熱段入口相對(duì)干度，hgl為汽化潛熱。

1.3方程的無(wú)量綱化

對(duì)上述方程中參數(shù)進(jìn)行無(wú)量綱化:

式中:L為長(zhǎng)度，cp為液相定壓比熱容，ΔTin為加熱段入口過(guò)冷度。將上述無(wú)量綱參數(shù)代入方程（1）～（3），得到如下的無(wú)量綱控制方程:

為了對(duì)上述的無(wú)量綱控制方程進(jìn)行簡(jiǎn)化，對(duì)其進(jìn)行量綱分析。結(jié)果發(fā)現(xiàn)由于Rgl數(shù)值接近0.001，導(dǎo)致αRgl等項(xiàng)與其余各項(xiàng)相比要小2～3個(gè)數(shù)量級(jí)，可以認(rèn)為是小量而忽略。此時(shí)控制方程簡(jiǎn)化為

同理對(duì)式（4）所示的輔助方程和式（5）～（8）所示的邊界方程無(wú)量綱化和量級(jí)分析，得到無(wú)量綱的輔助方程和邊界方程，具體過(guò)程不再冗述。

2　系統(tǒng)狀態(tài)空間表達(dá)

2.1線性增量方程的提出

假定系數(shù)具有線性特性，通過(guò)外加線性小擾動(dòng)的方法，對(duì)系統(tǒng)方程進(jìn)行線性化，得到描述系統(tǒng)動(dòng)態(tài)特性的線性增量方程。系統(tǒng)變量采用下式表示:

將式（16）代入式（13）～（15），得到線性增量的控制方程，其中單相區(qū)域的表達(dá)形式如下

而在兩相區(qū)域內(nèi)，線性擾動(dòng)方程的表達(dá)形式為

同時(shí)，將方程（16）代入無(wú)量綱的輔助方程，得到線性增量的輔助方程:

另外，將方程（16）代入無(wú)量綱的邊界方程中，得到線性增量的邊界方程:

2.2狀態(tài)空間方程組的建立

同時(shí)為了使變量和方程的數(shù)量相同，需要將線性增量的輔助方程式（24）代入控制方程，用以消去Δα這個(gè)變量，最終得到一組3×3的線性常微分方程組，即描述自然循環(huán)閃蒸不穩(wěn)定動(dòng)態(tài)特性的狀態(tài)空間表達(dá)式:

其中，矩陣M和N的參數(shù)如下所示，為了方便表述，首先定義如下的參數(shù):

此時(shí)矩陣M的參數(shù)可以表示為

而矩陣N的參數(shù)則可以表示為

2.3系統(tǒng)空間狀態(tài)方程穩(wěn)定域的判斷

對(duì)于方程（28）所示的狀態(tài)空間方程組，可以通過(guò)判斷系數(shù)矩陣N（=MTK）的特征根實(shí)部的符號(hào)來(lái)判斷系統(tǒng)的穩(wěn)定性。即求解如下的特征方程:

其特征根為

當(dāng)γi均小于零時(shí)，系統(tǒng)可以維持穩(wěn)定，即表示開(kāi)式自然循環(huán)系統(tǒng)能夠形成穩(wěn)定的閃蒸；而 γi存在大于零的情況時(shí)，系統(tǒng)無(wú)法維持穩(wěn)定，即此時(shí)的閃蒸流動(dòng)為閃蒸不穩(wěn)定工況。

3　模型計(jì)算結(jié)果及與實(shí)驗(yàn)結(jié)果比較

由于上述計(jì)算過(guò)程繁復(fù)，較難推導(dǎo)出閃蒸不穩(wěn)定區(qū)域所對(duì)應(yīng)的數(shù)學(xué)表達(dá)形式，因此本文利用MATLAB 7.0程序?qū)Ψ匠烫卣髦颠M(jìn)行求解和對(duì)系統(tǒng)穩(wěn)定性進(jìn)行判斷，其計(jì)算過(guò)程如圖2所示。

圖3給出了模型的計(jì)算結(jié)果，圖中實(shí)線表示計(jì)算得到的閃蒸不穩(wěn)定邊界，數(shù)據(jù)點(diǎn)表示開(kāi)式自然循環(huán)實(shí)驗(yàn)臺(tái)上得到的閃蒸不穩(wěn)定區(qū)域（實(shí)驗(yàn)裝置見(jiàn)圖4，實(shí)驗(yàn)方法見(jiàn)文獻(xiàn)［11］）。通過(guò)對(duì)比模型計(jì)算結(jié)果與實(shí)驗(yàn)結(jié)果發(fā)現(xiàn)，兩者之間趨勢(shì)上吻合較好。

其中在高加熱功率工況下，計(jì)算得到的不穩(wěn)定上邊界與實(shí)驗(yàn)結(jié)果存在一定偏差。此時(shí)實(shí)驗(yàn)觀察到加熱管出口過(guò)冷沸騰現(xiàn)象強(qiáng)烈而形成間歇噴泉流動(dòng)，而本模型中未考慮過(guò)冷沸騰現(xiàn)象對(duì)閃蒸流動(dòng)的影響，因此偏差產(chǎn)生的原因極可能是由于過(guò)冷沸騰過(guò)程的影響。

圖2　閃蒸不穩(wěn)定邊界的計(jì)算流程Fig.2 Calculation flow of flashing instability boundary

圖3　模型計(jì)算結(jié)果與實(shí)驗(yàn)結(jié)果的對(duì)比Fig.3 Comparison between the model calculation results and the experiment results

同樣在高加熱功率工況下，計(jì)算得到的不穩(wěn)定下邊界與實(shí)驗(yàn)結(jié)果也存在偏差。此時(shí)隨著加熱段入口過(guò)冷度降低，使得加熱段內(nèi)沸騰噴發(fā)過(guò)程增強(qiáng)，并與上升段內(nèi)閃蒸流動(dòng)過(guò)程共同影響自然循環(huán)的流動(dòng)特性［11］。而沸騰噴發(fā)過(guò)程會(huì)導(dǎo)致流量、上升段內(nèi)的溫度分布等出現(xiàn)劇烈地振蕩，進(jìn)而導(dǎo)致持續(xù)的穩(wěn)定閃蒸流動(dòng)無(wú)法形成，即可以認(rèn)為加熱段內(nèi)沸騰噴發(fā)過(guò)程對(duì)于穩(wěn)定閃蒸過(guò)程起到強(qiáng)烈的干擾作用。

因此本實(shí)驗(yàn)采用注氣誘導(dǎo)法排除沸騰噴發(fā)對(duì)于穩(wěn)定閃蒸過(guò)程的干擾作用。即在閃蒸不穩(wěn)定工況時(shí)，向上升段內(nèi)持續(xù)注入氣體，使得系統(tǒng)形成穩(wěn)定的空氣-水兩相流動(dòng)，此時(shí)自然循環(huán)流量快速增加，加熱段內(nèi)沸騰噴發(fā)過(guò)程被抑制。隨后逐漸減少注氣量直至停止，這一過(guò)程中閃蒸流動(dòng)過(guò)程會(huì)逐漸主導(dǎo)系統(tǒng)的流動(dòng)特性，而由于循環(huán)流量較高導(dǎo)致該過(guò)程中加熱段內(nèi)沸騰噴發(fā)不會(huì)出現(xiàn)。通過(guò)這種方法來(lái)排除沸騰噴發(fā)對(duì)閃蒸過(guò)程的影響，再通過(guò)觀察停止注氣后自然循環(huán)閃蒸流動(dòng)的發(fā)展演化過(guò)程來(lái)判斷閃蒸流動(dòng)是否穩(wěn)定。圖5給出了典型的注氣誘導(dǎo)法過(guò)程及其結(jié)果。

圖4 實(shí)驗(yàn)回路示意圖Fig.4 Sketch map of experimental loop

圖5（a）所示的實(shí)驗(yàn)工況位于模型計(jì)算的閃蒸不穩(wěn)定區(qū)域內(nèi)（圖3的不穩(wěn)定區(qū)域），注氣前系統(tǒng)處于閃蒸和沸騰共同誘發(fā)的流動(dòng)不穩(wěn)定狀態(tài)。當(dāng)注氣逐漸停止后，加熱管內(nèi)沸騰始終被抑制，而上升段內(nèi)閃蒸流動(dòng)出現(xiàn)振蕩，振幅不斷增大直至最終恢復(fù)為原來(lái)的流動(dòng)不穩(wěn)定狀態(tài)，即說(shuō)明此工況下閃蒸流動(dòng)過(guò)程為發(fā)散過(guò)程。而圖5（b）所示實(shí)驗(yàn)工況則處于模型計(jì)算的穩(wěn)定閃蒸流動(dòng)區(qū)域（圖3的下部穩(wěn)定區(qū)域），注氣前系統(tǒng)與工況1相似，但注氣逐漸停止后，上升段內(nèi)閃蒸流動(dòng)首先出現(xiàn)較弱振蕩，隨后振蕩逐漸減弱并形成穩(wěn)定閃蒸流動(dòng)，此時(shí)穩(wěn)定閃蒸驅(qū)動(dòng)的流量很大，使得加熱段出口為單相流動(dòng)，即工況2時(shí)自然循環(huán)閃蒸流動(dòng)過(guò)程為自發(fā)收斂過(guò)程，也證實(shí)了圖3實(shí)驗(yàn)中未出現(xiàn)穩(wěn)定閃蒸流動(dòng)主要是由于沸騰噴發(fā)過(guò)程的干擾作用。

圖6給出了通過(guò)重復(fù)實(shí)驗(yàn)得到的注氣誘導(dǎo)法形成的穩(wěn)定閃蒸流動(dòng)邊界，本文認(rèn)為是準(zhǔn)確的閃蒸不穩(wěn)定的區(qū)域下邊界，由圖可知該邊界與計(jì)算得到的不穩(wěn)定下邊界符合很好，即證明了本文模型結(jié)果具有很好的準(zhǔn)確性。另外，可以認(rèn)為本模型計(jì)算得到下穩(wěn)定區(qū)域內(nèi)的工況，都可以通過(guò)注氣誘導(dǎo)法將閃蒸不穩(wěn)定流動(dòng)轉(zhuǎn)化為持續(xù)穩(wěn)定閃蒸流動(dòng)。

圖5　注氣誘導(dǎo)法誘發(fā)穩(wěn)定閃蒸流動(dòng)的過(guò)程Fig.5 Flow rate variations of the method of gas injection to induce the steady flashing flow

圖6　注氣誘導(dǎo)法穩(wěn)定閃蒸邊界與模型計(jì)算結(jié)果的對(duì)比Fig.6 Comparison between model calculation results and experimental instability boundary by method of gas injection

4　結(jié)論

1）本文運(yùn)用狀態(tài)空間法對(duì)開(kāi)式自然循環(huán)閃蒸不穩(wěn)定進(jìn)行研究，提出了描述該過(guò)程的線性均相流模型，并導(dǎo)出了描述系統(tǒng)穩(wěn)定性的狀態(tài)空間方程。2）通過(guò)求解狀態(tài)空間方程得到了自然循環(huán)閃蒸不穩(wěn)定的邊界，將計(jì)算結(jié)果與實(shí)驗(yàn)結(jié)果對(duì)比發(fā)現(xiàn)兩者趨勢(shì)上吻合較好。其中上邊界的偏差是由于模型未考慮過(guò)冷沸騰的影響而造成的；下邊界的偏差是由加熱段內(nèi)沸騰噴發(fā)對(duì)穩(wěn)定閃蒸的干擾作用而產(chǎn)生的。3）通過(guò)注氣誘導(dǎo)法排除了沸騰噴發(fā)干擾作用，得到了準(zhǔn)確的閃蒸不穩(wěn)定下邊界，該邊界與計(jì)算得到的不穩(wěn)定下邊界吻合很好。4）通過(guò)實(shí)驗(yàn)發(fā)現(xiàn)，計(jì)算得到的下穩(wěn)定區(qū)域內(nèi)的工況點(diǎn)，均可通過(guò)注氣誘導(dǎo)法將閃蒸不穩(wěn)定工況轉(zhuǎn)化為穩(wěn)定閃蒸流動(dòng)工況。

［1］譚思超，龐鳳閣，高璞珍.自然循環(huán)過(guò)冷沸騰流動(dòng)不穩(wěn)定

性實(shí)驗(yàn)研究［J］.核動(dòng)力工程，2006，27（1）:18-21，93. TAN Sichao，PENG Fengge，GAO Puzhen.Experimental research on subcooled boiling flow instability of natural circulation［J］.Nuclear power engineering，2006，27（1）:18-21，93.

［2］KYUNG I S，LEE S Y.Periodic flow excursion in an open two-phase natural circulation loop［J］.Nuclear engineering and design，1996，162（2/3）:233-244.

［3］YANG S K.Stability of flashing-driven natural circulation in a passive moderator cooling system for Canadian SCWR［J］. Nuclear engineering and design，2014，276:259-276.

［4］吳莘馨，吳少融，姜?jiǎng)僖?，?閃蒸現(xiàn)象及其對(duì)自然循環(huán)穩(wěn)定性的影響［R］.中國(guó)核科技報(bào)告，1991. WU Xinxin，WU Shaorong，JIANG Shengyao，et al.Flashing in riser and its effect on flow stability in natural circulation system［R］.China nuclear science and technology report，1991.

［5］FURUYA F，INADA F，VAN DER HAGEN T H J J.Flashing-induced density wave oscillations in a natural circulation BWR-mechanism of instability and stability map［J］.Nuclear engineering and design，2005，235（15）:1557-1569.

［6］SAJITH P P，VEDULA R P，MITRA S K.An Investigation of flashing-driven natural circulation［J］.International journal of green energy，2006，3（4）:369-379.

［7］GUO Xueqing，SUN Zhongning，WANG Jianjun，et al.Numerical simulation of the transient behaviors in an open natural circulation system with a large scale［J］.Annals of nuclear energy，2015，77:83-93.

［8］GUO Xueqing，SUN Zhongning，WANG Jianjun，et al.Simulating investigations on the start-up of the open natural circulation system［J］.Nuclear engineering and design，2015，289:35-48.

［9］INADA F，F(xiàn)URUYA M，TASUO A.Thermo-hydraulic instability of boiling natural circulation loop induced by flashing （analytical consideration）［J］.Nuclear engineering and design，2000，200（1/2）:187-199.

［10］ZUBER N，F(xiàn)INDLAY J A.Average volumetric concentration in two-phase flow systems［J］.Journal of heat transfer，1965，87（4）:453-468.

［11］HOU Xiaofan，SUN Zhongning，F(xiàn)AN Guangming，et al. Flow characteristics in an open two-phase natural circulation loop［C］//Proceedings of the 22nd International Conference on Nuclear Engineering（ICONE 22）.Prague，Czech Republic，2014.

本文引用格式：

侯曉凡，孫中寧.開(kāi)式自然循環(huán)閃蒸不穩(wěn)定的線性均相流模型［J］.哈爾濱工程大學(xué)學(xué)報(bào)，2016，37（7）:930-935.

HOU Xiaofan，SUN Zhongning.Linear homogeneous-flow model of flashing instability in open natural circulation system［J］.Journal of Harbin Engineering University，2016，37（7）:930-935.

Linear homogeneous-flow model of flashing instability in open natural circulation system

HOU Xiaofan，SUN Zhongning
（Fundamental Science on Nuclear Safety and Simulation Technology Laboratory，Harbin Engineering University，Harbin 150001，China）

With respect to judging the flashing instability boundary in an open natural circulation system，we use a state-space method and develop a linear homogeneous-flow model to describe its physical process.From the statespace analysis results，we deduce a state space expression for estimating the stability of the system.We then determine the flashing instability boundary by solving these equations.The calculation result shows good agreement with the experimental result，with only a slight deviation under the condition of high heating rate.This deviation is mainly due to the interference effect of a boiling eruption in the steady flashing flow process in high heating power conditions，and this interference makes it difficult to establish a stable flashing flow.Next，we use a gas-injection method on the up-leg to avoid this interference effect in the experiments，and obtain an accurate lower boundary for the flashing instability.This boundary agrees well with the calculated lower boundary，which verifies the accuracy of this model.

open natural circulation；flashing instability；linear homogeneous flow model；state-space equations；gas injection method；steady flashing

10.11990/jheu.201506093

TL334

A

1006-7043（2016）07-930-06

2015-06-30.網(wǎng)絡(luò)出版日期：2016-05-06.

國(guó)家高科技研究發(fā)展計(jì)劃項(xiàng)目（2012AA050906）.

侯曉凡（1988-），男，博士研究生；孫中寧（1963-），男，教授，博士生導(dǎo)師.

孫中寧，E-mail:sunzhongning@hrbeu.edu.cn.

網(wǎng)絡(luò)出版地址：http：//www.cnki.net/kcms/detail/23.1390.u.20160506.1548.002.html